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Abstract 

The  semiconductor  laser  diode  offers  a  unique  system  to  investigate  nonlinear 
dynamics  when  optical  feedback  is  applied.  Although  there  is  extensive  research  of 
laser  diodes  with  optical  feedback  from  normal  dielectric  mirrors,  very  little  has  been 
done  experimentally  to  analyze  the  effects  of  degenerate  phase  conjugate  feedback 
from  a  BaTiOa  crystal.  This  research  experimentally  investigated  the  dynamics  of 
a  single-mode  laser  diode  with  weak  phase  conjugate  feedback  using  both  the  self- 
pumped  and  double  phase  conjugate  geometries.  The  experimental  results  validated 
a  mathematical  model  which  numerically  evaluates  the  Lang-Kobayahsi  coupled  dif¬ 
ferential  equations.  The  model  simulated  the  nonlinear  behavior  of  a  laser  diode 
subject  to  phase  conjugate  feedback  and  produced  simulated  RF  and  Fabry-Perot 
spectra  of  the  laser  output.  Model  and  experimental  results  showed  frequency  lock¬ 
ing  between  the  relaxation  oscillation  and  external  cavity  frequencies  and  changes 
in  the  relaxation  oscillation  as  a  function  of  cavity  length.  Resonant  behavior  in 
the  feedback  strength  necessary  to  undamp  the  relaxation  oscillation  frequency  as  a 
function  of  the  cavity  length  was  also  captured  experimentally  and  numerically.  Val¬ 
idation  of  the  model  and  experimental  results  presented  in  this  research  significantly 
contribute  to  the  understanding  of  the  nonlinear  behavior  of  a  laser  diode  subject  to 
optical  feedback. 
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Experimental  Investigation  of  Nonlinear  Dynamics  in  Single  Mode 
Semiconductor  Laser  Diodes  with  Phase  Conjugate  Feedback 


I.  Introduction 

Injecting  feedback  into  a  stable  oscillator  is  well  known  for  causing  the  oscil¬ 
lator  to  behave  nonlinearly.  In  particular,  a  laser  diode’s  response  to  feedback  has 
been  shown  to  include  highly  dynamical  and  even  chaotic  behavior  in  its  output.  The 
process  by  which  this  nonlinear  behavior  evolves  in  the  laser  is  not  well  understood 
and  is  the  subject  of  this  research.  This  dissertation  examines  a  semiconductor  laser 
diode  entering  into  a  chaotic  state  brought  about  by  the  application  of  phase  conju¬ 
gate  feedback.  The  investigation  was  performed  both  numerically  through  computer 
modeling  and  experimentally  using  two  different  methods  of  phase  conjugate  feed¬ 
back;  self-pumped  and  double  phase  conjugate  geometries.  This  document  starts 
with  a  brief  discussion  of  why  the  Air  Force  is  interested  in  chaotic  laser  diodes 
with  phase  conjugate  feedback  and  continues  with  an  in-depth  description  of  the 
theory  of  phase  conjugate  feedback  along  with  the  analysis  of  nonlinear  behavior 
used  throughout  the  research.  This  is  followed  by  a  discussion  of  the  computer  mod¬ 
els  used  in  the  research  and  experimental  descriptions  of  the  two  different  types  of 
phase  conjugate  geometries  investigated.  The  final  section  reviews  the  results  of  the 
experimental  and  model  outputs  and  draws  conclusions  from  the  research. 

1.1  Overview 

In  typical  optical  systems  using  laser  diodes,  the  laser  is  optically  isolated  from 
the  rest  of  the  system  to  prevent  feedback  from  unwanted  optical  reflections.  Optical 
feedback,  even  minute  amounts  in  semiconductor  laser  diodes,  on  the  order  of  -40  dB 
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can  easily  cause  the  laser  to  enter  into  a  chaotic  state  which  alters  the  single  mode, 
single  frequency  properties  of  the  laser,  thus  affecting  the  precise  control  of  the  laser 
in  typical  applications  such  as  coherent  communications  and  compact  disc  players. 

In  contrast,  the  Air  Force  and  other  agencies  are  investigating  methods  of 
exploiting  the  output  of  a  laser  diode  intentionally  placed  in  a  chaotic  state  (14)  (21). 
It  is  envisioned  that  the  broadband  and  nearly  random,  but  deterministic,  output, 
characteristic  of  a  laser  in  chaos,  will  provide  the  framework  for  a  new  generation 
of  secure  communications  for  air-to-air  and  ground-to-ground  networks.  It  is  this 
novel  application  that  forms  the  underlying  motivation  for  this  research.  A  complete 
understanding  of  this  dynamical  process  inside  the  laser  diode  and  the  coupling  to  its 
external  cavity  providing  the  optical  feedback  must  be  established  in  order  to  fully 
utilize  the  chaotic  state  of  a  laser  diode.  Theoretical  modeling  and  experimental 
work  have  been  conducted  by  the  Air  Force  and  civilian  institutions  on  laser  diodes 
intentionally  placed  in  a  chaotic  state  by  optical  feedback  from  conventional  dielectric 
mirrors,  commonly  referred  to  as  normal  feedback  (15)  or  through  the  use  of  optical 
injection  (14).  Only  recently  has  work  been  done  to  understand  the  role  of  phase 
conjugate  feedback  on  chaos  (26)  (2).  Research  conducted  under  this  effort  provides 
vital  details  for  understanding  the  onset  of  chaos  caused  by  phase  conjugate  feedback 
from  the  self-pumped  and  double  phase  conjugate  geometries. 

1.2  Introduction  to  Chaos 

What  is  Chaos?  When  a  system  is  said  to  have  chaotic  behavior,  what  kind 
of  behavior  does  this  entail?  In  terms  of  a  system’s  response,  chaos  can  be  defined 
as  aperiodic  long-term  behavior  of  a  system  whose  deterministic  outcome  is  highly 
sensitive  to  initial  conditions  (19)  (13).  In  other  words,  small  perturbations  to  the 
initial  conditions  can  cause  dramatically  different  outcomes  in  a  chaotic  system. 
The  field  of  study  of  this  complex  behavior  is  called  nonlinear  dynamics  where  a 
nonlinear  system  is  one  whose  dynamical  variables  describing  the  properties  of  the 
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system  (position,  velocity,  acceleration,  population,  amplitude,  etc.)  evolve  in  a 
nonlinear  fashion. 

For  the  laser  cavity  and  specifically  for  the  laser  diode,  chaos  can  be  achieved 
by  applying  optical  feedback  into  the  laser  cavity.  The  interference  of  the  interior 
optical  fields  with  the  delayed  externally  injected  fields  causes  the  laser  output  to  go 
through  a  series  of  changes  which  are  dependent  on  the  amplitude  and  phase  of  the 
optical  feedback.  With  extremely  small  amounts  of  feedback  the  laser  can  overcome 
the  effects  of  the  feedback  and  its  output  remains  consistent  with  a  solitary  laser. 
However,  at  a  critical  level  of  feedback,  the  stable  laser  output  will  change  to  a 
periodic  pulsed  output  and  finally  to  an  aperiodic,  or  chaotic,  output.  This  migration 
in  the  output  is  easily  shown  in  a  bifurcation  diagram  which  plots  amplitude  or 
intensity  output  against  feedback  level. 

Figure  1  shows  a  typical  bifurcation  diagram  of  a  fictitious  system.  In  this 
simple  example,  the  output  could  be  the  normalized  power  output  of  the  laser  mea¬ 
sured  at  regular  time  intervals  as  the  feedback  is  increased.  Initially  at  low  feedback 
levels,  the  laser  output  is  stable  with  a  solitary  lasing  frequency  Ui.  At  this  point  the 
feedback  is  too  weak  to  affect  the  laser  and  so  the  output  and  the  system  remains 
stable,  area  A.  As  the  feedback  level  is  increased,  a  Fold  bifurcation  occurs  and  the 
laser  hops  to  another  external  cavity  mode.  The  laser  continues  to  remain  stable, 
however,  it  now  lases  at  a  slightly  different  frequency  uj2.  The  effect  of  increasing 
the  feedback  in  region  B  now  increases  the  solitary  laser  output  amplitude  shown  as 
the  increasing  solid  line  for  area  B.  Further  increases  in  feedback  eventually  causing 
a  Hopf  bifnrcation  to  occur  and  the  output  amplitude  begins  to  oscillate  between 
two  values  as  shown  in  area  C.  In  this  region  the  relaxation  oscillation  frequency 
of  the  system  is  undamped  and  the  output  is  pulsating  at  this  frequency.  Further 
increasing  the  feedback  causes  the  output  to  go  through  another  Hopf  bifurcation 
creating  two  additional  branches.  Additional  increases  in  the  feedback  canse  more 
bifurcations  as  the  period  of  oscillation  doubles  (commonly  called  period  doubling). 
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Period  Doubled 
Periodic  Output 
Output  \ 


CHAOS 


Figure  1.  A  typical  bifurcation  diagram  showing  the  sudden  abrupt  changes  in  the 
behavior  of  a  fictitious  system  as  the  feedback  is  increased.  The  shaded 
area  denotes  the  area  of  chaotic  behavior.  The  output  axis  can  be  any 
measurable  change  in  the  laser  diode  system  such  as  output  power,  bias 
voltage,  bias  current,  etc.  The  letters  correspond  to  spectrums  shown  in 
Figure  2 
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Finally,  after  an  infinite  number  of  period  doubling  bifurcations,  the  output  becomes 
chaotic,  denoted  by  the  shaded  region  in  area  D  of  Figure  1.  The  frequencies  cre¬ 
ated  by  the  period  doubling  bifurcations  can  easily  be  seen  in  the  system’s  output 
spectrum.  Figure  2  shows  how  the  spectrum  of  this  fictitious  system  changes  as  the 
feedback  is  increased. 


-Vr  -Vr/2  0  vJ2 


-Vr  -Vr/2  0  Vr/2  Vr 


Figure  2.  Frequency  spectrum  of  the  fictitious  system  shown  in  Figure  1.  (A)  Stable, 
solitary  output  at  iOi,  (B)  New  stable  lasing  frequency  ui2,  (C)  Periodic 
output  where  the  relaxation  frequency  has  undamped  and  period  doubling, 
(D)  Chaotic  output  after  an  infinite  number  of  period  doublings. 


In  this  particular  example  the  process  leading  to  chaos  is  described  as  a  period¬ 
doubling  route  because  the  period  of  oscillation  doubles  after  each  bifurcation.  The 
period  doubling  route  is  common  in  semiconductor  laser  diode  experiments  with 
optical  feedback  and  extremely  short  external  cavities  (<1  cm).  However,  other 
quasi-periodic  routes  have  also  been  seen  (18).  Quasi-periodic  behavior  describes 
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the  condition  when  a  second  frequency  emerges  from  the  Hopf  bifurcation  and  the 
ratio  of  this  new  frequency  to  the  previous  frequency  is  not  rational. 

To  illustrate  a  few  of  the  key  terms  associated  with  chaos,  consider  a  simple 
one-dimensional  system  dX/dt  =  f{X).  Fixed  Points  for  a  system  are  values  X  such 
that  f{Xo)  =  0.  To  examine  the  behavior  of  the  system,  we  apply  a  Taylor  series 
expansion  to  the  equation  of  motion  around  the  fixed  point  X  =  Xg 


f{X) 


fiXg)  +  {X 
+  ^(X 


2^ 

dX^ 


+ 


(1) 


where  the  derivatives  have  all  been  evaluated  at  the  fixed  point  X  =  Xg-  If  we  let 
X  =  X  —  Xg  and  by  definition  f{Xo)  =  0  and  keep  only  the  first  term  of  the  expansion 
(a  valid  assumption  of  a  well-behaved  linear  system),  Eq.  (1)  can  be  rewritten  as 


fiX) 

dx 

dt 


dfjXg) 

dX 


X 


X 


for  which  the  solution  is 

x{t)  —  x{0)e^^ 


(2) 

(3) 

(4) 


where 


dfjXg) 

dX 


(5) 


is  the  Lyapunov  exponent  or  characteristic  value  of  the  fixed  point.  It  is  the  Lyapunov 
exponent  that  determines  the  behavior  of  the  system.  For  this  system,  A  has  three 
possible  values  (13:79): 


1.  A  <  0;  The  fixed  point  is  a  node  and  attracts  nearby  trajectories  (stable). 

2.  A  >  0;  The  fixed  point  is  a  repe//or which  repels  nearby  trajectories  (unstable). 
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3.  A  =  0;  The  fixed  point  could  be  a  node,  repellor  or  saddle  point.  It  is  necessary 
to  look  at  both  the  first  and  second  derivatives  of  x(t)  at  the  fixed  point  to 
determine  which  case  is  applicable.  For  a  node,  the  second  derivative  changes 
from  positive  to  negative.  The  opposite  occurs  for  a  repellor  for  A  =  0.  Saddle 
points  require  the  second  derivative  of  x{t)  to  be  the  same  sign  on  both  sides 
of  the  fixed  point. 

In  two  dimensional  (and  higher)  systems,  i.e.,  systems  with  two  or  more  inde¬ 
pendent  variables,  the  eigenvalues  can  include  a  complex  conjugate  pair  of  eigenval¬ 
ues  Ai_2  =  i?  ±  fO,  where  both  R  and  Q  are  real.  When  A  is  complex,  Eq.  (4)  can 
be  rewritten  to  see  the  eflfect  of  the  real  and  imaginary  parts  of  the  eigenvalue 


xit)  = 

=  a;(0)e^*  cos(Ot) 


(6) 

(7) 


We  see  that  the  trajectory  oscillates  at  frequency  and  is  either  damped  or  amplified, 
depending  on  the  sign  of  R.  For  R>  0,  the  trajectory  is  repelled  and  is  unstable  and 
when  R  <  0,  the  trajectory  is  attracted  to  the  fixed  point  and  is  stable.  When  the 
real  part  of  the  eigenvalue  is  0,  the  system  becomes  purely  periodic  which  results  in 
a  trajectory  that  forms  closed  loops,  or  cycles,  around  the  fixed  point.  Figure  3  is  a 
phase  plot  of  a  simple  harmonic  damped  oscillator  created  by  plotting  x{t)  versus  the 
velocity,  dx/dt,  to  show  the  behavior  of  the  trajectories  when  A  is  complex.  Table  1 
summarizes  the  behavior  of  the  system  for  the  different  eigenvalue  possibilities. 


Table  1.  Summary  of  the  behavior  of  trajectories  due  to  different  eigenvalues. 


Eigenvalue,  A 

Behavior 

Real,  Negative 

Real,  Positive 

Complex,  Negative  real  part 
Complex,  Positive  real  part 
Complex,  Real  Part  =  0 

Stable  Node 

Unstable  Repellor 

Spiral,  attracting  stable  focus 
Spiral,  repelling  unstable  focus 
Limit  Cycle 
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Figure  3.  Phase  plots  of  trajectories  of  a  simple  harmonic  damped  system  which 
has  a  complex  conjugate  pair  of  Lyapunov  exponents.  The  fixed  point 
for  this  system  is  shown  by  a  ‘+’  and  the  direction  of  increasing  time 
is  shown  by  the  arrow,  (a)  periodic  behavior  for  a  real  part  =  0,  (b) 
attracting  behavior  towards  the  fixed  point  when  the  real  part  <  0,  and 
(c)  the  trajectory  repels  away  from  the  fixed  point  out  to  the  limit  cycle 
or  spirals  into  the  limit  cycle  when  the  real  part  >  0 


According  to  the  Poincare-Bendixson  Theorem,  there  are  only  two  possible 
outcomes  for  the  trajectory  of  a  two-dimensional  system  with  complex  eigenvalues 
(13:104): 

•  The  trajectory  approaches  a  fixed  point  of  the  system  as  t  — >  oo.  (Figure  3b) 

•  The  trajectory  approaches  a  limit  cycle  as  t  ^  oo.  (Figure  3a,  c) 

The  impact  of  this  theorem  is  that  in  order  to  have  chaotic  behavior,  the  system 
must  be  of  dimension  three  or  higher  since  the  outcome  of  two-dimensional  behavior 
is  determined  by  the  Poincare-Bendixson  theorem  (13). 

An  excellent  example  of  a  two  dimensional  system  is  the  solitary  Class  B  laser 
(one  without  feedback).  This  system  can  be  modeled  as  a  simple  two-level  system, 
as  shown  in  Figure  4,  with  the  usual  pumping  rate,  Rp,  stimulated  transition  coef¬ 
ficients,  B21  and  Bi2,  and  the  spontaneous  transition  rate,  A21.  The  populations  of 
levels  1  and  2  are  Ni  and  N2,  respectively,  and  P  is  the  number  of  photons.  The  set 
of  coupled  differential  equations  describing  the  population  difference  between  levels 

-2BANP  +  2A2i{N  -  AN)  (8) 

BP  AN  -  7P  (9) 

where  B12  ~  P21  =  B,  N  =  {Rp/A2\  —  1)  Ni  is  the  unsaturated  population  of  level  2 
without  stimulated  emission  [P  =  0)  and  7  is  a  net  loss  of  photons  in  the  cavity  due 
to  various  loss  mechanisms  such  as  absorption  and  scattering. 

At  the  steady  state  condition  we  find  that  there  are  two  fixed  points: 

#1:  AN  =  N  ,  F  =  0 
#2:  AN  =  ^/B  ,  P  =  A2i(iV/7-l/P) 


1  and  2,  AN=  N2  —  Ni  are 

dAN 

dt 

d£^ 

dt 
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and  the  associated  Jacobian  matrix  for  the  eigenvalues  is 


^  -2^21  -  BP 
\ 


-2BAN 
BAN --y 


(10) 


When  the  Jacobian  is  evaluated  at  fixed  point  #1,  two  eignenvalues  can  be  found 


Ai  —  — 2^21 

A2  =  BN  —  7 


Since  A21  is  never  negative,  the  first  eigenvalue  is  always  negative.  Because  we  have 
no  emission  from  the  laser  at  P  =  0,  the  losses,  7,  are  greater  than  the  stimulated 
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emission.  Thus  the  second  eigenvalue  is  also  always  negative  for  this  fixed  point. 
Since  both  eigenvalues  are  negative  and  real,  fixed  point  7^1  is  a  stable  node  and 
describes  the  trivial  case  of  the  laser  operating  below  threshold.  However,  as  the 
laser  approaches  the  threshold  condition,  A2  becomes  less  negative  and  equals  zero 
precisely  at  threshold  since  at  this  point  the  gain  eqnals  the  losses.  At  this  precise 
condition,  fixed  point  #1  actually  becomes  unstable  as  the  laser  crosses  the  threshold 
condition.  As  we  shall  see,  this  unstable  condition  is  connteracted  by  fixed  point  #2 
which  becomes  stable  as  the  threshold  condition  is  reached  and  the  laser  operates  in 
a  continuous-wave  (c-w)  mode,  as  will  be  shown  below. 

Evalnating  the  Jacobian  at  the  second  fixed  point,  we  find  a  complex  conjugate 
pair  of  eigenvalues 


A±  = 


-BN  ^ 
7 


B^N^ 

(BN  A 

72 

A21 

( 7  V 

(11) 


Taking  the  Nd-YAG  laser  system  as  an  example.  Table  2  lists  typical  values  for 
the  parameters.  From  the  table,  it  is  easily  shown  that  {2y)/{A2i)  {BN » 
so  the  eigenvalues  can  be  simplified  to 


A±  = 


-BN 

7 


±  i 


27  (BN 
A21  \  7 


(12) 


For  lasing,  the  gain  is  always  larger  than  the  losses,  BN  >  7,  so  the  term  under  the 
radical  is  positive.  Thus  for  this  case,  fixed  point  #2  is  a  stable  spiral  node  since 
the  two  eigenvalues  have  negative  real  parts  and  non-zero  imaginary  parts.  Figure  5 
shows  a  graphical  view  of  this  fixed  point. 


Table  2.  Examples  of  typical  parameters  for  a  Nd-YAG  laser  system. 
A21  =  10"^  sec~^ 

N  =  7.33  *  10^8 
7  =  0.99  cm”^ 

B  =  2.1  —  8.8  *  cm^  /  sec 
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Population  Difference 


As  seen  in  this  figure,  the  damped  periodic  behavior  towards  the  steady  state 
output  value  is  the  relaxation  oscillations.  It  is  clear  from  this  simple  example  that 
a  system  with  only  two  degrees  of  freedom  cannot  reach  chaos  since  the  relaxation 
oscillations  are  always  damped  and  no  mechanism  exists  to  cause  them  to  become 
undamped.  Therefore  a  third  degree  of  freedom  is  necessary,  such  as  optical  feed¬ 
back,  for  chaotic  behavior  to  develop  as  stated  previously  by  the  PoincarABendixson 
Theorem. 

1.3  Phase  Conjugating  Mirrors  (PCM) 

To  provide  the  required  third  degree  of  freedom,  feedback  can  be  applied  to  a 
laser  diode  by  several  different  methods  including  electrical  modulation  of  the  bias 
current  and  injecting  optical  feedback  into  the  laser  cavity.  In  particular,  for  this 
research  a  phase  conjugating  mirror  is  used  in  order  to  examine  the  unique  effects 
of  the  external  cavity  formed  by  the  laser  diode  and  the  feedback  mirror.  PCMs 
are  unique  in  that  they  reflect  (produce)  a  phase  conjugate  beam  which  retraces 
the  exact  path  and  phase  front  of  the  input  beam.  Distortions  and  changes  induced 
by  optics  in  the  path  are  “un-done”  by  the  conjugate  beam  as  demonstrated  by 
Figure  6.  Therefore  at  any  point  along  the  optical  path  the  phase  front  of  the 
conjugate  beam  is  identical  to  that  of  the  original  beam.  The  two  types  of  PCMs 
used  in  this  research  are  created  by  a  BaTiOa  crystal  that  is  either  self  pumped  or 
double  pumped.  Figure  7  shows  a  schematic  view  of  the  two  configurations  used. 

The  self-pumped  CAT  geometry  shown  in  Figure  7  is  the  easiest  configuration 
for  establishing  a  conjugate  signal;  however,  it  has  the  disadvantage  of  having  the 
reflectivity  depend  solely  on  the  intensity  of  the  input  beam.  Si.  So  for  low  feedback 
levels,  the  input  beam  must  also  be  reduced  resulting  in  the  reduction  in  the  strength 
of  the  holographic  grating  formed  in  the  crystal.  But  this  slight  disadvantage  is 
greatly  out-weighed  by  the  fact  that  this  configuration  is  easy  to  implement  in  the 
lab  and  has  a  large  range  of  incident  angles  that  can  produce  an  adequate  conjugate 
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Figure  6.  The  reflected  beam  of  a  phase  conjugate  mirror  retraces  its  original  path, 
self-correcting  any  distortions  caused  by  optics.  The  lower  trace  shows  a 
beam  reflected  from  a  normal  dielectric  mirror  which  further  distorts  the 
wave-front  after  being  reflected. 
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Input  Beam 


Input  Beam 


Phase  Conjugate  of  B 


Phase  Conjugate  of  A 


Figure  7.  Two  different  methods  of  establishing  a  phase  conjugating  mirror  used  in 
this  research,  (a)  Self-pumped  CAT  geometry,  (b)  Double-phase  conju¬ 
gating  mirror  (DPCM). 


signal.  Garret  et  al  showed  that  a  range  of  incident  angles,  with  respect  to  the  c-axis, 
of  15°  -  75°  can  produce  reflectivities  as  high  as  40%  in  BaTiOs  at  a  wavelength  of 
514  nm  with  the  optimum  incident  angle  of  29.1°  (8).  The  CAT  geometry  also 
suffers  from  a  relatively  slow  rise  time  for  the  conjugate  signal  to  reach  90%  of  its 
steady  state  value.  Hengst  measured  this  rise  time  to  vary  between  <0.5  second  to 
tens  of  seconds,  depending  on  the  incident  intensity  and  angle  of  beam  Si  (12:22). 
However,  the  onset  time  to  establish  the  mirror  is  not  critical  for  this  research.  An 
investigation  of  the  dynamic  response  time  once  the  mirror  is  established  and  its 
effect  on  the  laser  dynamics  forms  a  part  of  the  present  investigation. 

The  transmission  grating  double  phase  conjugate  mirror  (DPCM)  configuration 
in  Figure  7b  offers  the  advantage  of  having  two  phase  conjugate  signals  to  use  as  a 
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feedback  signal  to  the  laser  diode  and  has  a  faster  onset  time  than  the  CAT  geometry. 
When  only  one  conjugate  signal  is  required  for  feedback,  the  second  beam  can  be 
blocked  with  a  high  quality  optical  isolator  placed  in  one  of  the  input  legs.  Weiss  et 
al  showed  that  in  BaTiOs  the  ratio  of  the  conjugate  intensities  is  equal  to  the  ratio 
of  input  intensities  (28) 

lAjout)  ^  iBjin) 

Isiout)  lA{in) 

HIa  (out)  is  chosen  as  the  feedback  signal,  this  configuration  also  allows  control 
of  the  conjugate  return  without  adjusting  the  input  intensity  /yi(m).  Now  the  re¬ 
flectivity  of  the  conjugate  mirror,  Rpcm-,  defined  as  the  ratio  of  the  output  intensity 
to  the  input  intensity,  is  proportional  to  the  second  input  beam,  RpcM  oc  Ipiin). 
This  allows  adjustment  of  the  reflectivity  without  affecting  the  input  beam  /^(w). 
In  addition,  this  conjugate  mirror  may  also  be  established  using  a  reflection  grating 
DPCM  configuration.  However,  current  research  shows  that  the  reflection  grating 
configuration  is  harder  to  implement,  and  may  not  offer  the  reliable  control  of  the 
conjugate  signal  required  for  the  experiments. 

Table  3  summarizes  the  dependencies  oi  Rpcm  and  the  controlling  parameters 
used  in  the  experiments.  For  each  configuration,  the  resultant  “reflected”  E-field  can 
be  described  as  a  modified  conjugate  of  the  input  field 

E*(t)  =  \l RpcME*{t)  exp{f(a;t  -I-  (14) 

where  E*{t)  is  the  complex  conjugate  amplitude  of  the  input  beam,  to  is  the  angular 
optical  frequency,  and  (t)pcM  is  a  phase  shift  of  the  conjugate  beam  relative  to  the 
probe  beam.  For  photorefractive  materials,  such  as  BaTiOs,  fbe  conjugate  beam  is 
in-phase  with  the  probe  beam,  therefore  (ppcM  —  0. 

As  a  starting  point  for  developing  the  theory,  the  next  chapter  review  the 
derivation  of  the  rate  equations  describing  the  laser  diode  system  without  feedback. 
The  phase  conjugate  term  discussed  earlier  in  this  section  is  then  applied  to  the 
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Table  3.  Listing  of  the  two  different  phase  conjugate  mirrors  and  their  controlling 
parameters  used  in  the  research. 


Mirror  Configuration 

Controlling  Parameter 

Dis/ Advantage 

CAT 

RpcM  ^  Si 

Reflectivity  values  High  but 
dependent  on  input  beam. 

Easy  alignment. 

DPCM 

RpcM  oc  lB{in) 

Reflectivity  values  Low  but 

can  be  varied  with  constant  signal 

beam. 

solitary  rate  equations  to  establish  the  three  coupled,  time-delayed,  rate  equations 
necessary  for  chaotic  development.  The  theory  section  is  then  followed  by  a  brief 
discussion  of  the  two  computer  models  that  numerically  integrate  the  set  of  nonlinear 
equations.  The  models  are  then  used  as  a  tool  to  analyze  the  experimental  data.  The 
remainder  of  the  dissertation  encompasses  the  experimental  descriptions,  experimen¬ 
tal  results  and  finally,  comparisons  between  the  model  and  analytical  predictions  and 
the  experimental  results. 
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IL  Theory  of  Nonlinear  Dynamics  for  Phase  Conjugate  Feedback  in 

Laser  Diodes 

Semiconductor  laser  diodes  are  extremely  sensitive  to  external  feedback  due  to 
the  low  reflectivities  of  their  output  mirrors  (facets).  The  low  reflectivity  allows  even 
minute  amounts  of  feedback  into  the  gain  region  of  the  laser  diode  which  can  cause 
the  diode  to  enter  into  a  state  of  chaos.  Optical  feedback  can  be  established  by  two 
methods,  conventional  feedback  (CF)  and  phase  conjugate  feedback  (PCF).  There 
has  been  extensive  research  investigating  conventional  feedback  with  normal  dielec¬ 
tric  mirrors  both  analytically  (22)(23)  and  experimentally  (16)(30).  In  this  research, 
phase  conjugate  feedback  is  used  for  the  study  of  nonlinear  dynamics  for  its  relative 
ease  of  implementation  and  its  unique  characteristic  phase  and  frequency  locking 
of  the  feedback  signal  not  found  in  other  types  of  feedback.  In  contrast  to  PCF, 
CF  requires  precise  alignment  and  focusing  of  the  feedback  beam  and  vibrations  in 
any  of  the  optics  can  induce  random  phase  variations  in  the  feedback  beam  which 
randomly  excite  various  external  cavity  modes.  Published  literature  on  the  study 
of  laser  diodes  with  PCF  has  focused  on  either  computer  models  that  show  chaotic 
behavior  (11)  or  experimental  observations  of  chaotic  behavior  with  PCF  (17).  How¬ 
ever,  there  is  a  definite  lack  of  published  material  that  validates  computer  models, 
and  analytical  approximations,  with  experimental  data.  Another  unique  aspect  for 
using  phase  conjugate  feedback,  is  that  there  is  no  conjugate  feedback  counterpart 
in  electrical  or  mechanical  feedback  systems  previously  studied  in  chaos. 

The  theory  of  phase  conjugate  feedback  begins  with  a  description  of  the  rate 
equations  for  the  dynamics  inside  the  laser  cavity  followed  by  the  inclusion  of  the 
phase  conjugate  feedback  term.  Figure  8  shows  a  diagram  of  the  different  components 
developed  in  this  section.  The  laser  cavity  on  the  left  emits  photons  into  the  external 
cavity  which  are  then  partially  retro-reflected  back  towards  the  laser  cavity  from  the 
phase  conjugate  mirror.  These  reflected  photons  reach  the  exit  facet  of  the  laser  with 


18 


19 


a  portion  re-entering  the  laser  cavity  and  the  rest  being  retro-reflected  back  into  the 
external  cavity  for  another  round  trip.  This  interface  at  the  laser  and  external  cavity 
sets  the  boundary  conditions  for  how  the  E-fields  interact  and  how  the  gain  medium 
is  altered  in  the  laser  cavity.  The  change  in  optical  gain  is  the  mechanism  that  drives 
the  laser  into  chaotic  development.  Once  a  complete  description  of  the  fields  inside 
the  laser  cavity  is  developed,  a  linear  stability  analysis  is  performed  to  determine 
operating  regions  and  feedback  levels  where  the  laser  diode  may  go  chaotic. 

2. 1  Development  of  Rate  Equations  Without  Feedback 

The  two  rate  equations  that  are  of  primary  concern  for  laser  diodes  are  the 
electron-hole  pair  number,  N(t),  and  the  optical  field,  E{t).  Starting  with  N{t),  Ya- 
mada  developed  an  equation  which  adequately  depicts  the  dynamics  of  the  electron- 
hole  pairs  in  the  gain  region  (29): 

f  =  (15) 

where  N{t)  is  the  number  of  electron-hole  pairs  and  DV'^N  accounts  for  the  diffusion 
of  carriers  in  the  gain  medium  (D  is  the  diffusion  coefficient).  Diffusion  of  carriers  can 
be  considered  instantaneous  (<C  1  ps)  on  the  time  scale  of  the  carrier  lifetimes  which 
are  typically  ~  1  ns,  and  therefore  can  be  neglected  in  the  case  of  semiconductor 
laser  diodes.  The  rate  equation  also  includes  the  addition  of  carriers  in  the  active 
region  of  length  d  due  to  an  applied  pumping  current  density,  J.  The  last  term, 
TZ{N),  is  the  lump  sum  rate  (#/sec)  of  all  loss  mechanisms  of  electron-hole  pairs 
which  include  both  stimulated  and  spontaneous  emission, 

N 

n{N)  =  Rst{N,P)  +  —  (16) 

tn 

Here  tn  is  the  electron-hole  carrier  lifetime  and  accounts  for  all  spontaneous  emission 
and  non-radiative  recombinations.  The  stimulated  radiation  rate,  Rst{N,P),  is  a 
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function  of  both  the  electron-hole  pair  number,  N,  and  the  laser  cavity  photon 
number,  P.  The  photon  number  is  a  dimensionless  quantity  and  represents  the 
magnitude  of  the  normalized  electric  field  with  the  energy  of  a  cavity  photon,  hui, 
such  that  P  =  \E\^.  Furthermore,  the  loss  of  a  carrier  pair,  5N  —  —1,  leads  to 
the  creation  of  a  photon,  6P  =  +1.  Van  Tartwijk  et  al  showed  that  the  stimulated 
emission  rate  can  be  written  as  (27) 


Rst{N,P)  =  v,g{N)P  (17) 

where  Vg  is  the  group  velocity  of  the  longitudinal  mode  inside  the  laser  cavity  and 
g{N)  is  the  electron-hole  dependent  gain  coefficient.  Based  on  numerical  calculations 
one  can  expand  the  gain  coefficient  around  the  threshold  value,  keeping  only  the 
lowest  two  terms  such  that 


g(N)  =  g(N,^)  +  ^(N  -  N„)  + 


(18) 


where  g{Nth)  is  the  threshold  gain  coefficient  and  Nth  is  the  number  of  electron-hole 
pairs  at  threshold  (27:101).  Since  the  electron-hole  pairs  are  restricted  to  remain  in 
the  active  region  and  the  injected  optical  field  may  extend  into  the  adjacent  areas 
outside  the  active  region,  the  resultant  optical  gain  of  the  laser  cavity,  gopt^  is  smaller 
than  g{N) 

9opt  =  r^f  (19) 


The  confinement  factor  or  filling  factor,  F,  describes  the  fraction  of  photons  that 
interact  in  the  gain  region  of  the  laser,  0  <  F  <  1.  Upon  substitution  of  Eqs.  (18) 
and  (19)  into  Eq.  (17),  the  stimulated  emission  rate  can  be  rewritten  as. 


Rst  — 


Vgg{Nth)  +  -  ^th) 


P 


^gQopt^Nth) 

f 


E^N-Nth) 


P 


(20) 
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r 


T+aN-  N,^) 


p 


where  To  =  Vggopt{Nth)  is  the  total  loss  rate  of  photons  in  the  laser  cavity  (s“^)  or 
commonly  called  the  inverse  photon  lifetime,  and  ^  =  Vgdg/dN  is  the  differential 
gain  rate.  Substituting  TZ{N)  and  Rst  into  Eq.  (16)  yields  the  expression  for  the  rate 
equation  for  the  electron-hole  pairs  inside  the  cavity  without  feedback, 


dN{t)  _  J  _  N{t) 
dt  e  Tjv 


^  +  {(W(f)  -  ]v,o 


m)f 


(21) 


where  E{t)  is  the  amplitude  of  the  complex  electric  field  E(f)  =  E{t)  exp^iut). 
Since  we  are  only  concerned  about  the  number  of  carriers  above  threshold,  we  can 
substitute  AN{t)  =  N{t)  -  Nth  and  let  p  =  (Jle)f(NthlTN)  be  the  ratio  of  carriers 
created  by  the  pumping  current  to  carriers  lost.  This  transforms  the  electron-hole 
rate  equation  into 


dANjt)^  Nth  ANjt) 

dt  ^  Thj  Tat 


+  iAN{t) 


(22) 


Development  of  a  similar  rate  equation  for  the  E-field  in  the  cavity  involves 
characterizing  the  gain  and  loss  mechanisms  for  the  active  layer  in  the  laser  diode 
and  examining  the  frequency  dependent  refractive  index.  Table  4  lists  the  primary 
loss  mechanisms  in  the  laser  cavity  that  are  considered  relevant  in  the  literature. 
The  material  loss  coefficient,  amati  is  derived  from  the  imaginary  part  of  the  laser 
gain  medium  susceptibility,  x,  for  the  unpumped  material  and  represents  photons 
lost  to  absorption  into  the  material  without  a  pumping  current  applied.  The  second 
loss  coefficient,  ainti  is  loss  of  photons  due  to  scattering,  with  a  being  the  scattering 
cross  section.  The  total  absorption  loss,  aabs,  is  the  sum  of  the  first  two  coefficients 
plus  an  additional  loss  of  photons  due  to  absorption  from  the  pumped  gain  medium. 
For  convenience,  the  optical  gain  coefficient  of  the  laser  diode  cavity,  gopt,  is  the  net 
gain  coefficient  in  the  semiconductor  gain  medium  which  incorporates  the  effect  of 
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absorption  from  the  unpumped  material  and  absorption  due  to  pumping: 

9opt  —  +  Xp)  (23) 

The  mirror  loss  coefficient  in  Table  4  is  the  effect  due  to  the  laser  facet  reflectivities, 
Ri  and  R2,  distributed  over  a  complete  round  trip  of  the  internal  laser  cavity,  2L. 


Table  4.  The  absorption  coefficients  (cm“^)  in  an  externally  pumped  FP-resonator 
consisting  of  an  active  medium  confined  between  mirrors  (27:98).  y"  and 
Xp  are  the  two  elements  of  the  total  scalar  susceptibility  in  an  isotropic 
material  due  to  unpumped  background  material  losses  and  losses  due  to 
external  pumping,  respectively. 


Material  Loss 

Oimat  = 

Internal  Loss 

^int  —  (J  /  TIC^.q 

Absorption  Loss 

^abs  —  (^mat  “h  ^int  T  ^"Kp/ 

Mirror  Loss 

am  =  (l/2L)ln(l/i?ii?2) 

To  derive  the  E-field  rate  equation,  we  start  with  the  round  trip  amplitude 
gain  threshold  condition  for  lasing 


Q  =  ^i?ii?2exp  |-2f^L  -  aabsL^  =  1 


(24) 


where  Q  is  the  round  trip  amplitude  gain  in  the  laser  cavity,  and  n  =n{N,  u)  is  the 
carrier  and  frequency  dependent  refractive  index  of  the  laser  gain  medium.  Since 
{Lon{N ,  (jj)) / c  is  carrier  and  frequency  dependent,  a  linear  expansion  around  the 
threshold  gives 


ton  ^  TT'th^th  ^ 
c  c 


coth  dn 
c  dN 


1  dn 

{N  —  Nth)  + 


^th) 


(25) 


where  Nth,  nth,  and  ujth  are  the  electron-hole  number,  refractive  index  and  optical 
frequency  at  the  laser  threshold,  respectively.  Substituting  Eq.  (25)  into  Eq.  (24) 
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yields 


g  =  |-2i  +  ‘^^{N  -  No.)  +  -  w,j)^  L  - 

(26) 

where  Ug  =  nth  +  cothdn/du)  is  the  group  refractive  index  for  the  longitudinal  mode 
of  the  laser  cavity.  Using  the  mirror  loss  term  from  Table  4,  this  equation  can  be 
re-written  in  exponential  form 

S  =  exp{-amL}  (27) 

(  2^  /  dfi  \  1 

X  exp  I  ---  inthUith  +  LOth-^iN  -  Nth)  +  -  uoth)  1  L  \ 

X  expj-of^ft^L} 


The  absorption  coefficient,  aabs^  can  also  be  written  as 

Olabs  —  <^int  +  ~(^o  +  Xp)  (28) 

—  ^int  9opt 

=  (Hint  -  Tg{N) 


where  Eq.  (19)  was  used  in  the  last  step.  Finally,  the  round  trip  amplitude  gain  is 
separated  into  the  product  of  a  frequency  independent  term,  Qi,  and  a  frequency 
dependent  term,  such  that  Q  =  QiQ^^  where 

=  expi^(^{Tg{N)-a,nt-arn)L-2i‘^^{N-Nth)'^'^  (29) 

g^  =  exp  -  2i^^{uj  -  a;tfe)|  (30) 

At  threshold,  N  =  Nth,  w  =  toth,  and  ^  =  1,  so  the  term  2nthUJthL/c  in  g^  must  be 
an  integer  multiple  of  27r  and  therefore  can  be  ignored.  Also,  using  Tin  =  2Lng/c  as 
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the  cavity  round  trip  time,  can  be  written  as 


^  ex.-p{-i(uj  -  uth)rin} 

=  exp{2n;tftri„}exp{-ia;ri„}  (31) 

Recognizing  that  iuj  in  the  Fourier  domain  for  a  monochromatic  field  yields  d/ dt  in 
the  time  domain,  this  may  be  rewritten  as 

Gu;  =  exp{iLOthTin}  exp{-Tind/ dt}  (32) 

Two  valid  assumptions  can  be  made  at  this  time.  The  first  is  that  the  frequency 

of  oscillation  is  close  to  the  threshold  value,  u  =  coth,  and  the  second  is  that  the  E- 

field  has  a  slowly  varying  amplitude  for  each  round  trip  inside  the  laser  cavity: 

E(t)  =  E{t)  exp{iutht}  (33) 

The  field  after  one  round  trip  in  the  gain  medium  is  then 

E(t  +  Tin)  =  GlGu,E{t)  =  Gi  exi>{iLOthTin}  exp  {-Tindidt}  E(t)  (34) 

The  derivative  in  the  exponential  is  expanded  as  a  power  series, 

exY>{-Tind/dt]E(t)  = 


g  {-Tind/dty 


n=0 

oo 


ni 


(35) 


E  ^  [{*  -  rin)  -  t]' 


n=0 


n\ 


E(t  -  Tin) 


which  applies  a  time  shift  to  the  E-field  such  that 

E(t)  =  g^E{t  -  Tin)  (36) 
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From  the  slowly  varying  amplitude  approximation  made  earlier,  the  time  shifted 
change  in  the  field  is  small  and  can  be  approximated  by  a  first  order  differential; 


E(t  -  Tin)  =  E(f)  -  Ti 


dEjt) 

dt 


(37) 


By  substituting  the  differential  approximation  in  Eq.  (37)  into  Eq.  (36),  the  E-field 
time  rate  of  change  now  becomes 


dE{t) 

dt 


^in 


l)m 


(38) 


Since  the  gain,  Gi,  is  close  to  unity  during  lasing,  an  exponential  series  expansion 
can  be  applied  to  the  term  (1  —  1/Gi)  such  that 


1  -  l/6l  «  (XaW  -  a,„,  -  a„)L  -  ^(/v  -  Na) 


(39) 


This  result  is  substituted  into  Eq.  (38)  producing  the  E-field  rate  equation 


dE{t)  _ 

dt  Tin 


Tg(N)  -  Cm  -  a„)  L  -  2iL‘^~AN(t) 


m 


(40) 


Substituting  for  the  full  expression  for  =  c/{2Lng)  and  Eq.  (18)  for  g{N)  gives 


E(^) 

dt 


{rg(N,i,)  -  o.:„,  -  a„)v,  +  Tv,^AN(t)  - 


E(t)  (41) 


At  the  threshold  lasing  condition,  the  cavity  gain  is  equal  to  the  sum  of  the  losses  in 
the  cavity,  and  so  the  first  term  in  brackets  vanishes.  Recalling  that  the  differential 
gain  rate  is  ^  =  Vgdg/dN,  the  E-field  rate  equation  can  be  simply  written  as 


dE{t) 

dt 


re- 


2iuJth  dn 


Tig  dN 


l\N{t)E{t) 


(42) 
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Finally,  the  substitution  of  a  as  the  line-width  enhancement  factor 


a  = 


cF 


dN 


(43) 


as  the  ratio  of  the  real  and  imaginary  parts  of  the  pump-induced  susceptibility,  Xp 
and  Xp  respectively,  the  E-field  rate  equation  can  be  written  as: 


(44) 


Equations  (22)  and  (44)  represent  the  set  of  coupled  rate  equations  for  the  semicon¬ 
ductor  laser  diode  without  feedback.  It  should  be  noted  that  these  equations  are 
only  valid  if  the  E-field  undergoes  small  variations  within  a  round-trip  time  (  2  ps) 

of  the  laser  cavity.  In  the  next  section  the  necessary  feedback  terms  are  added  to 
account  for  the  delayed  phase  conjugate  feedback  from  an  externally  pumped  phase 
conjugating  mirror. 


2.2  Rate  Equations  with  Phase  Conjugate  Feedback 

The  phase  conjugate  term  stated  earlier  can  be  used  directly  with  a  slight 
modification  to  account  for  the  time  lag  between  the  E-field  leaving  the  laser  cavity, 
and  re-entering  the  cavity  after  being  reflected  from  the  PCM, 

r  =  (45) 

c 

where  Lext  is  the  length  of  the  external  cavity  formed  by  the  laser  and  the  PCM. 
Adding  the  phase  conjugate  term  to  the  rate  equations  gives  the  full  coupled  rate 
equations  for  the  E-field  and  electron-hole  number  inside  the  laser  diode  cavity  cou¬ 
pled  to  an  external  phase  conjugating  cavity: 

=  ^(1  +  ia)Ti^N{t)E{t)  +  KW{t  -  t)  exp  [i(j>PCM  +  -  0  }  (46) 
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and  r]c  is  a  coupling  loss  term.  The  frequency  change  in  the  conjugate  beam  due 
to  nondegenerate  pump  beams  (up)  in  a  nondegenerate  four  wave  mixing  phase 
conjugate  mirror  is  Au)  =  u)p  —  u>.  Since  the  degenerate  case  is  used  in  this  research, 
Au)  =  0,  but  is  kept  for  generality. 

As  stated  earlier  chaotic  behavior  requires  systems  to  have  dimension  three  or 
higher.  It  is  convenient  to  split  the  complex  E-field  rate  equation  into  independent 
terms  for  the  normalized  amplitude  P{t)  =  |E(f)|  and  phase  The  final  set  of 
coupled  rate  equations  with  feedback  are: 


=  T(AN{t)P(t)+2K\/P(t)P(t  -  t)  cos  (t>{t  -  r)  -h  4){t)  -\-  (pPCM  +  2Aui{t  - 

(49) 

=  iar^AiV(t)  -  +  <^(0  +  4>pcm  +  2Auj{t  - 

(50) 


dAN{t) 


=  (P-  1)- 


AN{t) 


-  ^+eAiv(^)]  iP(t)p 


2.3  Regions  of  Laser  Operation 

The  roadmap  for  investigating  the  nonlinear  dynamics  of  the  coupled  rate  equa¬ 
tions  involves  finding  the  steady  state  solutions,  performing  a  linear  stability  analysis 
and  examining  the  solutions  as  the  feedback  level  is  increased.  The  eigenvalues  for 
the  system  can  either  be  real  or  a  complex  conjugate  pair.  A  real  eigenvalue  that  is 
identically  =  0  (Fold  bifurcation)  implies  that  the  overall  phase  of  the  laser  output  is 
indeterminate  which  is  typical  for  feedback  with  a  normal  dielectric  mirror  (normal 
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feedback).  In  normal  feedback  cases,  the  phase  of  the  laser  is  strongly  dependent 
on  the  placement  of  the  feedback  mirror  and  small  variations  in  the  external  cavity 
length  can  cause  the  phase  to  vary  wildly.  When  a  fold  bifurcation  occurs  the  laser 
is  actually  mode  hopping  from  one  eigenvalue  to  another  giving  the  appearance  that 
the  laser  is  exiting  chaos  to  a  stable  state.  For  PCF,  a  real  eigenvalue  of  A  =  0 
does  not  occur.  It  will  be  shown  later  that  a  laser  with  PCF  will  phase  lock  and  be 
completely  deterministic  even  for  small  amounts  of  feedback.  For  the  situation  of  a 
complex  conjugate  pair  of  eigenvalues,  a  Hopf  bifurcation  occurs  when  the  real  part 
is  equal  to  zero  as  it  switches  from  negative  (damped)  to  positive  (undamped).  The 
imaginary  part  of  the  eigenvalue  gives  the  new  frequency  of  oscillation  undamping 
at  this  instability  which  amplitude  modulates  the  lasing  frequency. 

Table  5.  With  increasing  feedback,  the  laser  output  can  be  categorized  into  five 
different  levels  of  chaotic  development. 

Level  I  Initial  reduction  of  laser  linewidth 

Level  II  Relaxation  oscillations  become  undamped 

Level  III  Excitation  of  external  cavity  modes 

Level  IV  Significant  broadening  of  laser  linewidth.  Chaos 

Level  V  Linewidth  narrows  and  remains  narrow  for  further  increases  in  feedback 


The  main  focus  of  this  research  is  the  dynamics  of  the  laser  leading  up  to,  but 
not  entering,  chaotic  activity.  It  has  been  widely  accepted  that  a  chaotic  system 
can  be  categorized  into  five  characteristic  levels  of  chaotic  development.  The  first 
of  these  levels.  Level  I,  is  identified  by  an  initial  reduction  in  the  laser  linewidth, 
which  typically  occurs  with  feedback  levels  <  0.1%  of  the  laser  output.  As  the 
optical  feedback  in  increased.  Level  II  is  reached  when  relaxation  oscillations  become 
undamped  and  begin  to  appear  in  the  frequency  spectrum.  Further  increases  in  the 
feedback  causes  the  laser  to  enter  Level  III  and  are  marked  by  the  excitation  of 
external  cavity  modes.  At  Level  IV,  sometimes  called  coherence  collapse  or  chaos, 
the  laser  linewidth  broadens  significantly  from  a  few  MHz  to  several  GHz.  Finally, 
Level  V  is  reached  when  the  laser  exits  out  of  chaos  with  a  reduced  linewidth,  similar 
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to  Level  I.  Level  V  is  different  from  Level  I  in  that  it  needs  extremely  strong  feedback 
levels  and  usually  requires  the  exit  facet  of  the  diode  to  be  anti-reflection  (AR) 
coated  to  increase  the  coupling  of  the  feedback  beam  into  the  laser  cavity.  Table  5 
summarizes  the  five  levels.  In  this  research  Levels  I,  II,  and  III  are  investigated  in 
detail. 
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III.  Mathematical  Analysis 


The  three  coupled  time-delayed  differential  equations  presented  at  the  end  of 
the  previous  chapter  are  extremely  difficult  to  analyze  analytically.  This  chapter 
discusses  the  mathematical  techniques  and  approximations  that  are  used  to  find 
solutions  to  the  equations.  The  analytical  approximations  made  in  this  chapter 
will  be  revisited  in  the  discussion  of  the  experimental  results  in  order  to  check  the 
validity  of  the  approximations.  In  addition  to  the  analytical  tools  presented  here, 
two  computer  models  that  numerically  evaluate  the  three  coupled  equations  are 
introduced  in  the  next  chapter. 

3.1  Scaling 

To  perform  numerical  analysis  and  computer  modeling  on  the  three  coupled 
nonlinear  Eqs.  (49)-(51),  it  is  beneficial  to  carry  out  a  scaling  to  non-dimensionalize 
the  equations.  The  scaling  used  in  this  research  is  the  same  as  Gavrielides  et  al  (9) 
which  gives: 

Time  normalization  to  the  photon  life  time:  s=  t/jp 
Feedback  normalization:  r]=  KTp 
Carrier  decay  time  normalization:  T=  r^/Tp 
Field  normalization:  A  = 

Carrier  Normalization:  N=  ^T^ANrp 

The  resulting  equations  with  the  new  normalized  parameters  are: 

^  =  aN  +  r]a(s  —  r)  cos  [0(s  —  r)  +  (l)(s)  +  (t)pcM  +  2Auj{s  —  r/2)]  (52) 

as 

^  =  aN  -  sin  [0(s  -  r)  +  0(s)  +  (t)pcM  +  2Aa;(s  -  t/2)]  (53) 

ds  a[s) 


31 


(54) 


T-—  =  V-N  -{I  +  2N)a^ 
as 

where  ^4(5)  =  r  has  been  normalized  with  Tp  and  V  =  (p  —  l)Nth  is  the 

amount  of  pumping  above  threshold.  For  the  experiments,  the  totally  degenerate 
condition  Aoj  =  0  is  used  but  Au  is  retained  here  for  generality. 

3.2  Steady  State  Solutions 

A  common  technique  for  examining  the  behavior  of  the  system  is  to  look  in 
the  steady  state  where  all  derivatives  can  be  set  to  zero  and  examine  the  solutions 
and  the  behavior  of  the  system.  Substituting  the  steady  state  values  for  a,  N  and  (p: 

a  =  as  N  =  Ns  (p  —  ut  +  (ps  (55) 

where  as,  Ns  and  (ps  are  all  constants  and  u)  is  the  angular  frequency  of  the  laser 
measured  relative  to  the  solitary  frequency  w,  yields  the  solution 

uj  —  Auj  (56) 

=  0  for  degenerate  PCF. 

This  simple  equation  indicates  that  the  relative  frequency  Co  is  identical  to  the  fre¬ 
quency  difference  between  the  pump  beam  and  the  solitary  laser.  In  other  words,  it 
shows  that  the  laser  locks  to  the  frequency  of  the  pump  beam  of  the  phase  conjuga- 
tor.  This  condition  offers  a  unique  ability  of  the  PCM  to  adjust  the  lasing  frequency 
through  adjustments  of  the  pump  frequency.  A  similar  adjustment  of  the  laser  fre¬ 
quency  with  normal  feedback  is  not  possible.  In  addition  to  the  frequency  locking 
condition  given  above,  we  also  have 

Cj  =  —rjy/l  -\-  sin('!?  +  arctan(Q;))  (57) 


32 


where  ■!?  =  20^  +  (f)pcM-  This  solution  imposes  phase  locking  of  the  laser.  Notice 
that  this  solution  is  independent  of  r  which  is  also  a  major  change  from  the  case  of 
the  normal  mirror.  For  a  normal  mirror  the  phase  is  a  strong  function  of  the  delay 
time  T  and  therefore  minor  changes  in  the  external  cavity  length  can  adversely  affect 
the  phase.  By  using  phase  conjugation,  the  r  dependence  is  eliminated  and  Eq.  (57) 
requires  the  phase  of  the  laser  to  lock  to  a  constant  value  set  by  the  conjugate  mirror 
4>pcMi  the  material  parameter  a,  and  the  feedback  strength  rj.  From  this  equation, 
it  is  clear  that  the  solution  for  only  exists  for  ]a)|  <  {ri\/l  +  ct^),  where  only 
the  positive  slope  is  stable.  This  equality  describes  the  fold  bifurcation  line  in  the 
(r},u})  plane.  Figure  9  shows  the  graphical  solution  to  Eq.  (57)  and  the  saddle-node 
bifurcation  lines  dividing  the  space  into  three  regions.  Points  below  the  straight 
lines  fall  into  the  region  of  periodic  solutions  which  have  been  associated  with  four- 
wave  mixing  (FWM)  solutions.  Solutions  above  the  bifurcation  lines  are  the  phase 
and  frequency  locked  solutions  given  by  Eq.  (57).  The  relevant  parameters  and  a 
MathCad^  template  for  the  figure  can  be  found  in  Appendix  A. 


Taking  into  account  the  restrictions  given  by  the  inequality,  the  solutions  for 


-d  are: 

■d  =  —  arctan(Q;)  —  arcsin 
for  the  nodes,  and 


U) 


^??vtT 


-f  2mr  n  =  0, 1,  2, . . . 


(58) 


7?  =  —  arctan(Q;)  —  arcsin  (  — .  )  -f  {2n  -1-  l)7r  n  =  0, 1,  2, . . .  (59) 

V77V1  + 

for  the  saddles.  These  solutions  are  also  shown  graphically  in  Figure  9a. 

In  addition,  the  steady  state  electron-hole  number  and  photon  number  are 


Ns  =  —r]  cos'd  (60) 


^MathCad  is  a  registered  trademark  of  MathSoft 
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1 

2 


(61) 


f  V  +  r]  cos  7?  \ 

Y 1  —  2?/  cos 'd  j 

The  numerator  in  Eq.  (61)  represents  the  net  rate  of  electron-holes  between  pnmping, 
V,  and  the  steady  state  population  difference  Ng-  The  denominator  is  related  to  the 
stimulated  emission  gain  given  in  Eq.  (21)  for  emission  of  photons  due  to  stimulated 
emission.  The  effect  of  the  feedback  is  to  change  the  threshold  gain  rate  of  the  laser, 
G=  1  —  2?7COS'd,  with  the  sign  of  cosi9  being  the  dominant  factor.  In  particular,  the 
threshold  gain  is  reduced  for  I'd]  <  7r/2  with  the  strongest  reduction  occurring  when 
=  0  or  Au  =  —rja.  This  condition  is  unique  to  non-degenerate  phase  conjugation 
where  the  difference  between  the  pump  and  probe  frequencies  can  be  adjnsted  to 
meet  the  Au  =  —r]a  requirement.  For  example,  given  typical  values  of  rj  —  0.012 
and  o;  =  4,  a  separation  of  the  pump  and  probe  frequencies  by  1.836  GHz  would 
produce  the  lowest  threshold  gain,  G  =  1—27].  Thus,  operating  at  this  point,  the  laser 
with  a  non-degenerate  phase  conjugate  mirror  would  extract  the  optimum  gain  from 
lasing  medium  since  the  gain  threshold  is  reduced  to  the  minimum  level.  Of  course, 
this  also  implies  that  for  other  phase  conjugate  configurations  where  Auj  ^  —7]a,  and 
particularly  for  the  degenerate  case,  Aco  =  0,  the  phase  of  the  laser  will  not  operate 
at  the  point  of  lowest  threshold  gain,  but  rather  at  a  point  yielding  the  smallest  line 
width  (1)(17)(20). 

3.3  Steady  State  Stability 

It  was  stated  earlier  that  the  steady  state  solutions  to  the  coupled  nonlinear 
equations  were  independent  of  r  for  the  phase  conjugate  case.  Although  the  steady 
state  solutions  may  be  independent  of  the  external  cavity,  the  stability  of  these  so¬ 
lutions  could  depend  on  r.  To  investigate  the  sensitivity  of  the  steady  state  solution 
to  small  perturbations,  a  linear  stability  analysis  can  be  performed  on  the  rate  equa¬ 
tions  by  linearizing  around  their  steady  state  values  and  assuming  a  time  dependence 
of  the  form  exp(5r).  The  result  is  a  secular  determinant  equation  whose  solutions. 
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]j(^s)  =  0,  yields  the  time  dependent  behavior  of  the  perturbation: 

D(s)  =  [Ts  +  {1 +  2al)][s^ +  2r]s cos'd  +  7]“^ {1  - 

+  2('P  +  ?]cos'd)[s  +  r7\/T+~^cos('d  +  arctan(Q;))(l  +  e~®'^)]  (62) 

From  the  fold  bifurcation  lines  in  Figure  9,  it  is  easily  identified  that  the  fold  bifurca¬ 
tion  for  the  totally  degenerate  case  occurs  for  zero  feedback,  rj  —  0.  Since  this  point 
offers  no  additional  information,  it  is  the  Hopf  bifurcation  point  that  is  of  interest 
for  the  dynamics  of  the  laser  for  degenerate  phase  conjngate  feedback.  At  the  Hopf 
point,  the  real  part  of  the  complex  pair  of  roots  of  Eq.  (62)  passes  from  negative  to 
positive  values  and  is  set  to  zero.  Substituting  s  =  iCl  into  Eq.  (62)  and  separating 
the  real  and  imaginary  parts  gives  two  equations  for  the  feedback  level  r],  and  the 
frequency  fl  of  the  resulting  periodic  solution 

-Q,[2'r]fl  cos  d  +  tf  sin2Qr]  -I- e- - - - ^[^?^(1  -cos2f2r)  -  + 

^  1  -  277  cos  w 

2e{V  +  7]  cos  d)ri\/l  +  q;^cos(7?  -|-  arctan(Q;))(l  -I-  cos  Or)  =  0  (63) 


and 


e - [27?f2  cos  d  +  7]^  sin  2fir]  -|-  f2[77^(l  —  cos  2D,r)  —  -f 

1  —  277  cos  d 

2e{V  4-  77  cos  7?)  [n  —  'q'Jl  +  cos{d  -f  arctan(Q;))  sin  fir]  =0  (64) 

where  e  =  l/T.  At  this  point  these  two  equations  must  be  solved  with  a  nonlinear 
root  finder  to  find  exact  solutions,  or  approximations  must  be  made  so  that  analytical 
solutions  can  be  found. 

Nonlinear  behavior  of  systems  at  the  Hopf  bifurcation  point  provides  insight 
into  the  magnitude  of  and  77.  Typically  the  frequency  emerging  out  of  the  Hopf 
bifurcation  is  close  to  the  relaxation  frequency  of  the  system  or  one  of  its  harmonics. 
In  addition,  from  the  behavior  of  Eqs.  (63)  and  (64)  the  combination  of  Vtr  must  be 
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0(1).  Furthermore,  rj  is  typically  small  and  is  0(e).  This  leads  to  an  appropriate 
scaling  of 

D.  =  e^^^a  T  =  r]  =  eE  (65) 

where  a,  0,  and  E  are  all  0(1).  Substitution  into  Eq.  (63)  and  (64)  leads  to  a 
solution  for  the  feedback  and  frequency  at  the  Hopf  point  in  the  limit  as  e  ^  0  (9) 


_ e(l  +  2V) _ 

2[\/r+~^cos('d  +  arctan(£r))  cos^(\/^©/2)  —  cos'd] 


+  0(e='/2) 


n  ^  c  ^  ^  arctan(ft))  sin(v^0)  ^  cjc^'^) 

^  4[\/l  +  cos('d  +  arctan(Q;))  cos2(\/^0/2)  —  cos'd] 

(67) 

where  the  lowest  order  in  e  is  the  most  dominant  term  and  is  the  relaxation 
frequency  that  emerges  from  the  Hopf  bifurcation. 


Gavrielides  pointed  out  that  this  perturbation  expansion  leading  to  Eqs.  (66) 
and  (67)  fails  when  cos(V^^0/2)  approaches  0  which  occurs  at  the  resonant  points 
of  Q.T  =  (2n  +  l)7r  or 

—  =  +  n  =  0,l,2,...  (68) 

^ext  ^  / 

Equation  (68)  is  the  condition  for  frequency  locking  between  the  relaxation  oscillation 
and  external  cavity  frequency  uiext  with  the  locking  ratio  being  (n+|).  Equations  (67) 
and  (68)  also  show  that  when  the  Hopf  frequency  equals  the  solitary  relaxation 
frequency  of  the  laser,  Qh  =  lor,  the  external  cavity  length  (r)  falls  in  the  middle  of 
the  locking  region  for  adjacent  n. 


An  additional  simplification  that  can  be  used  is  the  case  of  a  ^  oo.  This 
condition  has  been  extensively  examined  and  is  well  understood  with  externally 
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injected  laser  systems  and  offers  a  more  compact  form  for  rjn  and  D.h' 


1  +  2V  u)R 
2V  2cos2(^) 


(69) 


and 

=ujr(i-  ^  tan  )  (70) 

where  ujr  =  y^(2P)/r  is  the  scaled  relaxation  frequency  of  the  solitary  laser.  These 
two  equations  provide  a  look  into  the  behavior  of  both  r]H  and  as  the  external 
cavity  length  is  varied.  Figure  10  shows  both  functions  over  a  range  of  = 
18-28  cm.  The  unbounded  behavior  of  the  functions  is  due  to  the  lack  of  a  loss 
mechanism  in  the  system  and  the  fact  that  these  equations  are  only  valid  in  the  non¬ 
resonance  lengths.  The  measured  resonant  peaks  from  the  experimental  data  will  be 
bounded  by  inherent  loss  mechanisms  and  most  likely  will  have  lower  magnitudes. 
Nonetheless,  the  general  behavior  provides  an  excellent  preview  of  the  expected 
behavior  for  r)H  and  VLh- 

The  analytical  techniques  and  approximations  presented  in  this  chapter  offer 
one  method  of  analyzing  the  set  of  three  nonlinear  equations.  The  next  chapter  dis¬ 
cusses  the  computer  models  which  numerically  evaluate  the  equations  to  hnd  steady 
state  solutions.  The  computer  models  have  been  developed  to  simulate  the  develop¬ 
ment  of  chaotic  behavior  of  the  laser  diode  and  are  used  interpret  the  experimental 
data. 
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IV.  Computer  Modeling  of  the  Phase  Conjugate  Feedback 

The  previous  chapter  provided  insight  into  the  necessity  of  using  computer 
modeling  to  examine  the  nonlinear  behavior  of  time-delayed  coupled  differential 
equations.  Computer  modeling  provides  an  alternative  to  analytical  solutions  with¬ 
out  making  approximations,  other  than  the  basic  assumptions  leading  to  Eqs.  (52)- 
(54).  This  chapter  contains  a  brief  discussion  of  the  two  models  used  in  support  of 
this  research.  The  models  were  developed  by  the  Nonlinear  Optics  Center,  Phillips 
Laboratory,  Kirtland  Air  Force  Base,  New  Mexico. 

4.1  Bifurcation  Model 

The  LKCONBI2.C  model  integrates  the  Lang-Kobayashi  (L-K)  phase  conju¬ 
gate  equations  to  produce  a  listing  of  data  points  which  can  be  plotted  to  form 
bifurcation  diagrams.  A  copy  of  the  C-source  code  is  included  in  Appendix  B.  In¬ 
puts  to  the  model  include  external  cavity  length,  beginning  and  ending  values  of 
feedback,  laser  pumping  above  threshold,  linewidth  enhancement  parameter,  and 
lifetime  ratio  T.  The  scaling  used  in  the  model  is  the  same  as  presented  in  Section 
3.1. 

One  minor  change  to  the  original  LKCONBIF.C  program  was  made  to  move 
the  user  input  parameters  to  a  separate  input  file.  This  was  done  primarily  to  avoid 
compiling  the  source  code  for  multiple  changes  in  the  input  parameters  and  to  incor¬ 
porate  additional  non-scaled  experimental  parameters  such  as  measured  relaxation 
oscillation  frequency  and  external  cavity  length.  A  listing  of  all  the  input  parameters 
is  shown  in  Table  6. 

The  LKCONBI2  program  actually  contains  two  separate  models:  the  evalua¬ 
tion  of  the  full  L-K  rate  equations  presented  in  Fqs.  (52)-(54)  and  the  evaluation  of 
the  Phillips  Laboratory’s  phase  equation  (3).  The  user  may  select,  via  an  input  pa¬ 
rameter,  which  model  to  evaluate.  Both  operate  the  same  throughout  the  program. 
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Table  6.  Input  parameters  contained  in  a  separate  ASCii  file  for  the  LKCONBI2 
model. 


Variable  Name 

Purpose 

BEGIN 

Beginning  integration  time 

END 

Ending  integration  time 

EXTIME 

Time  when  steady-state  begins 

NPOINTS 

Number  of  point  to  use  in  integration 

KMAX 

Integration  resolution  magnifier 

HTAMIN 

Minimum  feedback  level 

HTAMAX 

Maximum  feedback  level 

HTASTEP 

Number  of  steps  between  max  and  min  feedback 

pi 

Pumping  level  above  threshold 

alpha 

Linewidth  enhancement  factor 

extJength 

Unsealed  external  cavity  length  (cm) 

delta 

Scaled  frequency  offset  from  solitary  laser 

phipcm 

Phase  shift  due  to  PCM 

photonJife 

Unsealed  photon  lifetime  (sec) 

carrierJife 

Unsealed  carrier  lifetime  (sec) 

model 

l=full-up  equations,  0=phase  equation 

Y 

Gain  saturation  term 

except  for  the  equation(s)  used  in  the  integration  routine.  For  simplicity,  the  full 
rate  equation  option  will  be  discussed  below. 

The  model  begins  by  dividing  the  full  range  of  feedback  specified  by  the 
user  into  feedback  steps,  (HTAMAX-HTAMIN)/HTASTEP.  Increasing  the  value 
of  HTASTEP  provides  more  detail  in  the  bifurcation  diagrams  but  also  increases  the 
run-time  of  the  model.  Typical  values  of  HTASTEP  are  50  to  100.  Small  values 
of  HTASTEP  are  not  recommended  since  bifurcation  points  could  be  missed  with 
too  large  a  step  size.  For  each  step  of  the  feedback  parameter,  a  time  series  of  the 
system  is  reconstructed  using  a  fourth-order  Runge-Kutta  integration  routine  on 
the  set  of  normalized  differential  equations,  Eqs.  (52)  -  (54).  The  delayed  values  of 
the  amplitude,  phase  and  carrier  number  needed  by  the  equations  are  handled  by 
establishing  a  circular  stack  in  memory,  whose  length  is  four  times  the  normalized 
external  cavity  round-trip  delay,  r.  An  index  into  the  stack  is  used  to  retrieve  the 
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proper  values  and  the  stack  is  initialized  with  zeros  so  that  the  feedback  does  not 
take  effect  until  at  least  one  round  trip  in  time  r  has  been  completed.  The  factor 
of  four  in  the  stack  size  is  to  account  for  the  fourth-order  integration  routine  which 
needs  a  delayed  value  for  each  integration  step. 

Extrema  of  the  amplitude  are  not  saved  until  the  user-established  steady  state 
time,  EXTIME,  has  elapsed.  This  allows  transients  in  the  time  series  to  diminish 
and  steady  state  conditions  recorded.  During  the  steady  state  condition,  when  a 
peak  in  the  amplitude  has  occurred,  the  feedback  rate,  time,  amplitude  and  phase 
of  the  laser  output  are  recorded  in  the  data  file.  These  points  can  then  be  plotted 
using  a  separate  plotting  program  to  create  bifurcation  diagrams.  Figure  11  shows 
two  examples  of  bifurcation  diagrams  plotting  the  output  amplitude  extrema  versus 
feedback  level.  The  relevant  parameters  for  these  plots  are  T  =  316.1,  o;  =  4,  and 
V  =  3.0. 

External  cavity  lengths  can  be  divided  into  two  regimes  based  on  the  value  of 
the  external  cavity  and  relaxation  frequencies.  An  external  cavity  is  labeled  as  being 
short  when  Vext  >  while  long  cavities  have  v^xt  <  i^r-  Each  external  cavity  regime 
has  separate  characteristics  which  are  also  found  into  the  computer  models.  Figure  11 
shows  a  computed  bifurcation  diagram  for  the  short  cavity  (a)  and  long  cavity  (b). 
The  main  differences  between  the  short  and  long  regimes  occurs  at  the  Hopf  points. 
In  the  short  cavity,  the  bifurcation  at  the  Hopf  points  is  distinct  indicating  that  only 
one  frequency  is  emerging  at  a  time.  This  distinction  is  due  in  part  to  the  fact  that 
the  external  cavity  frequency  is  greater  than  the  relaxation  oscillation  frequency  and 
requires  a  higher  rj  to  undamp.  The  frequencies  that  emerge  at  additional  Hopf  points 
are  related  to  the  first  Hopf  point  through  a  period  doubling  or  quasi-periodicity  of 
the  relaxation  frequency  that  has  already  undamped.  On  the  other  hand,  the  Hopf 
bifurcations  in  the  long  cavity,  beyond  the  first  bifurcation,  are  more  obscured.  A 
spectral  analysis  using  the  LKCON2.C  model  described  in  the  next  section  shows 
that  the  external  cavity  frequency  has  undamped  which  is  quickly  followed  by  several 
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Figure  11.  Bifurcation  diagrams  created  using  the  LKCONBIF.C  model,  (a)  Short 
external  cavity  with  r  =  16,  extJength  =  1.0  cm,  (b)  Long  external 
cavity  with  r  =  640,  extJength  =  40  cm. 
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harmonics.  This  process  then  develops  immediately  into  a  chaotic  state  shown  in 
the  bifurcation  diagram. 

A  direct  one-to-one  comparison  between  the  bifurcation  model  and  experimen¬ 
tal  data  was  not  intended  in  this  research.  To  produce  a  bifurcation  plot  experi¬ 
mentally,  a  data  acquisition  system  to  record  the  amplitude  time  series  of  the  laser 
output  would  require  sampling  the  laser  amplitude  at  a  minimum  of  twice  the  re¬ 
laxation  oscillation  frequency.  For  the  laser  at  hand,  vr  ~  3.7  GHz  at  low  power, 
so  a  sampling  frequency  of  8  GHz  would  be  required.  The  data  acquisition  system 
available  for  the  experiments  had  a  time  series  measurement  limit  of  500  MHz  and 
therefore,  most  of  the  experimental  data  consisted  of  RF  and  optical  spectra. 

Although  the  experimental  equipment  could  not  meet  the  fast  sampling  re¬ 
quirements  to  generate  bifurcation  diagrams,  the  LKCONBI2.C  model  was  still  use¬ 
ful  in  providing  insight  to  bifurcation  points  and  the  type  of  route  to  chaos  due  to 
various  changes  to  experimental  parameters.  The  most  effective  tool  for  comparing 
theoretical  and  experimental  data  was  the  spectrum  model  discussed  in  the  next 
section. 

4-2  Spectrum  Model 

The  second  computer  model  is  the  LKCON2.C  program  whose  C-source  code 
is  listed  in  Appendix  C.  This  model  is  used  to  create  RF  and  FM  spectra,  which 
can  be  compared  to  the  RF  spectrum  analyzer  and  Fabry-Perot  (F-P)  outputs.  The 
majority  of  the  model  is  similar  to  the  LKCONBI2.C  program  except  that  a  single 
value  for  the  feedback  parameter  is  used  and  the  number  of  points  used  in  the 
integration  routine,  NPOINTS,  is  increased  to  provide  a  larger  bandwidth  in  the 
FFT  algorithm.  The  original  Phillips  Laboratory  LKCON.C  program  was  modified 
in  a  similar  fashion  as  the  LKCONBIF.C  model  to  move  the  user  input  parameters 
to  a  separate  ASCii  file  to  simplify  multiple  runs  of  the  model  with  different  input 
parameters.  The  user  input  parameters  are  the  same  as  those  listed  in  Table  6  with 
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the  addition  of  an  extra  term,  inject,  which  is  the  feedback  level  where  the  spectra 
are  to  be  calculated. 

From  the  bifurcation  diagram,  the  user  can  choose  a  value  of  the  feedback  to 
use  in  the  LKCON2.C  model.  Once  the  integration  of  the  nonlinear  equations  is 
complete  for  the  feedback  level,  two  FFTs  are  applied  to  the  reconstructed  time  se¬ 
ries.  The  first  one  is  applied  to  the  amplitude  of  the  signal  and  gives  the  expected  RF 
spectrum.  The  second  FFT  is  applied  to  the  phase  term  and  gives  the  FM  spectrum. 
Both  are  output  in  separate  data  files,  FOURIERRF.DAT  and  FOURIERFM.DAT, 
respectively.  The  computed  time  series  data  is  also  recorded  in  the  LKCON.DAT  file. 
Figure  12  is  an  example  of  the  RF  and  FM  outputs  of  the  LKCON2.C  model.  The 
bifurcation  map  is  the  same  as  Figure  11b  and  the  spectra  have  been  calculated  for  a 
feedback  rate  of  rjinject  =  0.0065.  The  RF  spectrum  shown  in  Figure  12b  represents 
the  external  cavity  frequencies  and  higher  harmonics.  This  spectrum  also  shows 
peak  spacings  of  185  MHz  and  375  MHz  due  to  a  quasi-periodicity  of  the  external 
cavity  frequency.  Spectra  similar  to  Figure  12b  will  be  compared  with  spectra  from 
the  RF  spectrum  analyzer  during  the  experimental  portion  of  this  research.  The 
FM  spectrum  in  Figure  12c  shows  the  relaxation  oscillation  frequency  at  ±  4.7  GHz 
along  with  modulations  from  the  RF  surrounding  ur  and  the  DC  component.  This 
type  of  spectrum  is  representative  of  the  Fabry-Perot  output  in  the  experiments.  In 
a  latter  section,  the  model  outputs  and  experimental  spectra  are  compared  using  RF 
and  FM  spectra  for  numerous  external  cavity  lengths. 

As  with  the  bifurcation  model,  this  model  behaves  differently  in  the  short  and 
long  cavity  regimes.  For  short  cavities  where  >  ^Ri  the  spectra  are  consistent 
with  period  doubling  routes  to  chaos.  The  dominant  frequency  in  this  case  is  the 
relaxation  oscillation  which  quickly  undamps  and  period  doubles  with  increasing 
feedback.  The  external  cavity  frequency  has  little  effect  on  the  chaotic  development 
since  the  feedback  level  required  to  undamp  Vext  is  significantly  higher  than  the 
relaxation  oscillation.  At  these  large  levels  of  feedback,  the  relaxation  oscillation 
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has  already  period-doubled  many  times  and  the  laser  begins  to  operate  in  Level  IV, 
chaos.  Although  the  short  cavity  region  could  not  be  directly  tested  experimentally, 
the  conclusions  drawn  on  the  behavior  of  the  model  are  consistent  with  theory  and 
experiments  using  short  cavities  with  normal  dielectric  mirrors  and  direct  optical 
injection  (24). 

In  the  long  cavity  regime  where  Uext  <  ^r-,  the  model  as  described  earlier 
runs  into  a  difficulty  in  predicting  the  external  cavity  frequency  and  mode  spacing. 
The  result  is  an  external  cavity  frequency  that  ranges  from  half  the  experimental 
value  of  Vext  for  the  lower  end  of  the  long  cavities  to  one-fourth  v^xt  for  very  long 
cavities.  This  problem  was  investigated  and  found  to  be  caused  by  a  fundamental 
limitation  of  the  model  which  did  not  include  gain  saturation  of  the  laser.  The  model 
was  corrected  to  include  a  Gain  Saturation  term,  Y,  which  yields  results  that  more 
closely  agree  with  experimental  data.  The  specifics  of  the  changes  made  to  the  model 
are  discussed  in  Section  6.1. 

The  LKCON2.C  model  can  also  be  used  to  predict  the  steady  state  behavior 
of  the  laser  and  phase  conjugate  feedback  system  as  a  function  of  the  external  cavity 
length.  Figures  13-15  show  model  predictions  of  the  dependency  of  the  external 
cavity  frequency,  relaxation  oscillation  and  locking  ratio  as  functions  of  external 
cavity  length.  These  plots  will  be  compared  with  experimental  results  to  validate 
the  performance  of  the  model,  as  shown  in  Chapter  VI.  It  is  clear  from  these  figures 
that  the  model  output  is  representative  of  the  trends  introduced  by  the  analytic 
approximations  in  the  previous  chapter. 
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V.  Experimental  Descriptions  and  Results 

The  experimental  portion  of  this  research  can  be  divided  into  three  groups  of 
testing.  The  hrst  set  of  experiments  characterized  the  solitary  laser  without  feed¬ 
back  to  document  the  operating  characteristics  of  the  Spectra  Diode  Laser  (SDL) 
5412  laser  diode.  In  the  second  group  of  experiments,  weak  feedback  from  the  CAT 
phase  conjugate  mirror  geometry  was  applied  to  the  laser.  Measurements  of  the  laser 
output  due  to  the  effects  from  the  conjugate  feedback  were  conducted  by  varying 
the  feedback  power,  laser  bias,  and  external  cavity  length.  The  final  group  of  exper¬ 
iments  used  the  double  phase  conjugate  mirror  (DPCM)  geometry  as  the  feedback 
mechanism  and  performed  similar  tests  as  with  the  CAT  mirror  for  comparison. 
Each  group  of  experiments  is  discussed  in  detail  along  with  the  results.  This  chapter 
is  followed  by  an  analysis  of  the  comparison  between  model  and  experimental  data. 

5. 1  Solitary  Laser  Characterization 

Before  feedback  experiments  were  conducted,  it  was  necessary  to  measure  some 
key  characteristics  of  the  solitary  laser  to  see  the  effects  due  to  chaos.  The  parameters 
measured  in  these  simple  experiments  without  feedback  included  lasing  threshold, 
light-current  output,  nominal  wavelength,  and  linewidth.  The  equipment  and  con¬ 
figuration  used  in  this  test  are  shown  in  Figure  16. 

The  lasing  threshold  was  established  by  increasing  the  laser’s  bias  current  from 
zero  mA  and  collecting  output  power  measurements  just  after  the  microscope  ob¬ 
jective  (MO)  collimating  lens.  Power  measurements  were  taken  using  a  Newport 
815  power  meter  and  the  lasing  threshold  was  extrapolated  from  the  point  in  the 
light-current  curve  where  the  power  output  suddenly  increases.  After  threshold,  the 
bias  current  was  further  increased  in  steps  to  its  maximum  value  of  85  mA  so  that 
the  efficiency  or  slope  of  the  light-current  line  could  be  found.  Figure  17  shows  the 
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light-current  data  for  the  solitary  laser.  The  extrapolated  threshold  current  is  found 
to  be  14.92  mA  with  a  slope  efficiency  of  0.93  mW/mA. 

In  a  similar  experiment,  a  0.6  m  monochromator  and  optical  multichannel 
analyzer  (OMA)  combination  was  used  to  measure  the  center  wavelength  and  de¬ 
termine  if  the  laser  was  truly  single  mode  as  the  bias  current  was  varied  from  just 
past  threshold  to  approximately  three  times  threshold.  As  shown  in  Figure  16,  at 
the  same  time  the  OMA  was  measuring  the  center  wavelength,  a  F-P  was  used  to 
measure  the  linewidth  and  quantify  the  amount  of  noise  on  the  laser  waveform.  The 
main  source  of  noise  was  due  to  current  and  temperature  fluctuations  from  the  power 
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Table  7.  Summary  of  the  solitary  laser  c 

laracteristics. 

Lasing  threshold 

14.92  mA 

Light-Curve  Slope 

0.93  mW/mA 

Nominal  center  wavelength  (40  mA  bias) 

804.9  nm 

Change  in  wavelength  vs  bias 

0.03  nm/mA 

Nominal  Linewidth 

135  MHz 

supply  which  led  to  temporal  linewidth  fluctuations  in  the  F-P  display.  Figures  18- 
20  show  OMA  scans  of  the  solitary  laser  biased  at  30,  40  and  50  mA.  At  a  bias  of  30 
mA,  Figure  18,  the  OMA  shows  a  nominal  center  wavelength  of  804.4  nm  but  also 
some  weak  peaks  surrounding  the  main  lobe.  Although  this  small  structure  is  not 
contributing  a  significant  amount  of  power  to  the  solitary  laser  output,  it  was  found 
in  later  feedback  experiments  that  these  other  modes  become  undamped.  The  F-P 
spectrum  shown  in  Figure  21  with  the  laser  biased  at  30  mA  also  depicts  frequencies 
other  than  the  main  lobe  at  -1.5  and  -2.2  GHz.  Increasing  the  bias  current  to  40 
mA,  Figure  19,  reduces  second  order  modes  and  the  laser  operates  single  mode  with 
a  center  line  of  804.9  nm  and  a  nominal  linewidth  of  135  MHz.  At  a  bias  current  of 
50  mA,  Figure  20,  the  laser  is  still  single  mode  with  a  center  wavelength  805.1  nm. 
Because  of  the  uncertainty  of  single  mode  operation  for  laser  bias  currents  below 
40  mA,  the  feedback  experiments  described  in  the  next  several  sections  were  con¬ 
ducted  at  biases  of  40  and  45  mA,  or  2.6x  and  3x  threshold.  Table  7  summarizes 
the  characteristics  of  the  solitary  laser. 

5.2  Self-Pumped  Phase  Conjugate  Feedback 

Figure  22  shows  the  experimental  setup  for  investigating  the  laser  dynamics 
using  the  CAT  geometry  as  the  feedback  mirror.  The  major  components  of  this 
configuration  were  a  SDL  5412  single  mode  laser  diode  (LD)  and  the  self-pumped 
phase  conjugate  mirror  created  by  a  2400  ppm  Rh-doped  BaTiOs  crystal.  The  laser 
output  was  collimated  using  a  20X  (8  mm  focal  length)  microscope  objective  then 
immediately  split  with  a  50/50  BK-7  non-polarizing  beam  splitter  (BS).  The  reflected 
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Figure  20.  Optical  spectrum  of  the  solitary  laser  at  bias  current  50  mA. 
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Figure  21.  Fabry-Perot  output  of  the  solitary  laser  biased  at  30  mA.  The  linewidth  of  the 
laser  line  is  measured  to  be  130  MHz.  The  center  lobe  of  the  scan  is  rough  due 
to  low  frequency  noise  in  the  power  supply  and  temperature  instabilities  in  the 
T-E  cooler.  The  smaller  peaks  at  -1.5  GHz  and  -2.2  GHz  are  other  frequencies 
that  exist  at  low  bias  currents. 
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Figure  22.  Diagram  of  the  setup  used  in  the  CAT  experiments 
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light  from  the  beamsplitter  was  passed  through  an  anamorphic  prism  pair  (APP) 
to  re-shape  the  beam  into  a  circular  form  which  was  well  suited  for  the  apertures 
in  the  diagnostic  equipment.  An  optical  isolator  (01),  -41  dB,  was  used  to  prevent 
any  back  reflections  from  the  diagnostic  equipment  from  re-entering  the  laser  diode 
thus  establishing  an  unwanted  external  cavity.  To  examine  the  laser  output  for  RF 
frequencies,  an  Electro- Optics  Technology  (EOT)  amplified  p-n  junction  diode  (PIN 
diode),  frequency  response  >  800  MHz,  was  used  with  a  Hewlett-Packard  (H-P)  RF 
Spectrum  analyzer  which  had  a  frequency  cut-off  of  1.8  GHz.  Two  Fabry-Perots 
(F-P)  were  available  to  examine  frequencies  beyond  the  1.8  GHz  cut-off  of  the  RF 
spectrum  analyzer.  One  F-P  had  an  8  GHz  FSR  and  Finesse  >  200,  the  second 
with  a  6  GHz  FSR  and  Finesse  >  10000.  The  former  proved  to  be  more  useful  in 
examining  the  optical  spectrum  of  the  laser  diode  since  multiple  F-P  modes  could 
not  be  completely  eliminated  in  the  latter.  The  final  piece  of  diagnostic  equipment 
was  the  0.6  meter  monochromator/OMA  which  was  used  to  verify  that  single  mode 
operation  was  maintained  during  the  experiments. 

For  the  beam  going  to  the  crystal,  a  half-wave  plate  (HW)  was  used  to  rotate 
the  nearly  linearly  polarized  laser  output  to  coexist  with  the  extraordinary  direction 
in  the  crystal  which  is  the  favored  orientation  for  generating  a  phase  conjugate  signal. 
In  later  experiments,  the  half-wave  plate  was  eliminated  to  reduce  the  number  of 
optical  elements  in  the  external  cavity  so  that  shorter  external  cavity  lengths  could 
be  obtained.  This  did  not  have  an  adverse  effect  on  the  conjugate  mirror  since  the 
diode  output  is  greater  than  95%  linearly  polarized  and  the  laser  is  orientated  so 
its  primary  polarization  is  in  the  extraordinary  direction  for  the  BaTiOa  crystal. 
A  series  of  neutral  density  filters  and  a  variable  neutral  density  filter  were  used  to 
control  the  amount  of  feedback  and  a  50  mm  focal  length  lens  (LI)  was  used  to  focus 
the  beam  loosely  into  the  crystal.  The  phase  conjugate  return  beam  was  sampled 
with  a  Newport  815  power  meter  to  measure  the  average  power  returned  to  the  laser 
from  the  external  cavity  mirror  and  a  second  EOT  photodiode  could  also  be  used  in 
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place  of  the  Newport  meter  so  that  the  RF  spectrum  of  the  return  signal  could  be 
analyzed. 

5.2.1  Single  Mode  Operation  with  Phase  Conjugate  Feedback.  The  primary 
purpose  for  this  experiment  was  to  verify  that  the  laser  remained  in  single-mode 
operation  during  the  onset  of  chaos.  If  the  laser  did  not  remain  single-mode  then  the 
assumption  of  a  single-mode  laser  would  have  to  be  corrected  and  the  coupled  set  of 
nonlinear  Eqs.  (49)-(51)  would  have  to  be  modified  to  include  multi-modes.  As  was 
mentioned  in  an  earlier  section,  the  OMA  and  F-P  scans  showed  possible  multi-mode 
operation  of  the  solitary  laser  at  laser  current  bias  levels  of  30  mA.  Although  the 
multiple  modes  in  the  solitary  laser  spectrum  at  that  bias  level  did  not  contribute 
significant  power  to  the  laser  output,  applying  feedback  had  a  dramatic  affect  on 
these  small  amplitude  modes  causing  them  to  significantly  increase  in  amplitude. 
As  mentioned  earlier,  the  feedback  experiments  must  be  conducted  with  the  laser 
operating  in  a  single  mode  so  an  appropriate  bias  level  to  maintain  single  mode 
operation  was  required. 

The  OMA  was  used  to  monitor  the  laser  output  as  the  feedback  was  increased 
from  zero  /rW  until  the  laser  began  to  show  signs  of  the  relaxation  oscillation  and 
external  cavity  frequencies  undamping.  Figure  23  shows  outputs  of  the  OMA  with 
the  laser  biased  at  30  and  40  mA  with  weak  phase  conjugate  feedback.  For  30 
mA,  a  second  mode  appears  at  804.6  nm  for  extremely  weak  feedback,  <  0.002  % 
of  the  laser  output,  even  before  the  relaxation  oscillation  frequency  was  detected. 
Additional  peaks,  as  shown  in  the  figure,  appear  at  804.7  and  804.9  nm  for  weak 
to  moderate  feedback.  The  inset  in  the  figure  is  the  RF  spectrum  taken  at  the 
same  time  as  the  OMA  scan  which  shows  i'ext=  225  MHz  and  higher  harmonics 
have  become  undamped.  When  the  laser  is  biased  at  40  mA,  the  laser  remains 
single  mode  during  the  onset  to  chaos.  Figure  23b  shows  the  OMA  scan  when  the 
laser  is  biased  at  40  mA  with  moderate  feedback  (>  0.01%).  The  inset  shows  the 
strong  RF  peaks  corresponding  to  Uext  and  higher  harmonics.  Additional  increases 
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in  the  feedback  from  that  shown  in  the  figure  resulted  in  a  broadband  RF  structure 
(chaotic  behavior).  Therefore,  at  the  time  of  the  scans  shown  in  the  figure,  the  laser 
was  operating  just  below  the  chaotic  region  which  covers  the  full  range  of  feedback 
levels  investigated  in  this  research.  Hence,  laser  bias  currents  above  40  mA  were 
suitable  for  the  additional  feedback  experiments  and  ensured  single  mode  operation 
during  weak  to  moderate  feedback  levels. 

5.2.2  CAT  Phase  Conjugate  Mirror  Time  Response.  Part  of  the  charac¬ 
terization  of  the  onset  of  chaos  for  the  CAT  mirror  configuration  was  to  measure  the 
speed  of  the  crystal  as  a  phase  conjugate  mirror.  This  would  determine  if  a  fourth 
coupled  differential  equation  depicting  the  time  response  of  the  reflectivity  of  the 
PCM  was  necessary  to  describe  the  nonlinear  behavior  of  the  laser  diode  completely. 
The  time  dependent  behavior  of  the  mirror  reflectivity  would  manifest  itself  in  a 
dn/dt  equation.  The  theory  leading  up  to  the  set  of  coupled  nonlinear  equations, 
Eqs.  (49)-(51),  has  assumed  that  the  PCM  performs  almost  instantaneously,  so  k  is 
considered  a  constant. 

The  context  of  speed  of  the  PCM  is  different  than  the  onset,  or  build-up  time 
of  the  mirror.  The  onset  is  the  time  needed  by  the  crystal  to  create  a  conjugate 
signal  once  the  crystal  is  illuminated.  In  BaTiOs  and  other  photorefractives,  the 
holographic  grating  constituting  the  PCM  is  established  by  photoionization  of  im¬ 
purities  creating  charge  carriers.  These  charge  carriers  then  migrate  in  the  crystal 
material  and  are  trapped  in  the  dark  regions  of  the  light  intensity  pattern.  This 
sets  up  a  space-charge  electric  field  which  in  turn  modulates  the  refractive  index  of 
the  crystal  through  the  Pockels  effect.  Since  the  whole  process  of  establishing  the 
holographic  index  grating  relies  on  the  migration  of  charge  carriers  in  the  crystal,  the 
fundamental  time  limit  for  the  onset  of  the  grating  is  on  the  order  of  2  ms  for  BaTiOs 
(31:102).  However,  once  the  holographic  grating  has  been  established,  the  speed  at 
which  the  mirror  reacts  to  small  changes  in  the  input  beam  is  the  real  question  at 
hand.  The  development  of  chaotic  behavior  in  the  laser  diode  introduces  numer- 
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ous  high  frequencies  in  its  output  resulting  from  the  undamping  of  the  relaxation 
oscillations  (several  GHz)  and  external  cavity  frequencies  (several  MHz). 

These  questions  were  answered  in  a  series  of  experiments  that  compared  the 
frequency  content  and  delay  of  the  feedback  signal  from  the  CAT  phase  conjugate 
mirror  configuration  to  a  feedback  signal  from  a  normal  dielectric  mirror.  Figure  22 
shows  the  experimental  setup  for  the  CAT  phase  conjugate  mirror.  Fast  PIN  diodes 
were  used  to  sample  the  direct  laser  output  and  the  feedback  signals.  The  output 
signal  from  the  PIN  diodes  was  displayed  on  a  RF  signal  analyzer  and  a  500  MHz 
oscilloscope.  In  both  cases  the  feedback  was  adjusted  so  that  the  relaxation  oscilla¬ 
tion  and  the  external  cavity  frequencies  became  slightly  undamped  so  that  the  laser 
would  not  enter  the  chaotic  regime.  Furthermore,  the  intensities  of  the  beams  enter¬ 
ing  the  PIN  diodes  were  balanced  with  a  variable  neutral  density  (ND)  filter  placed 
in  front  of  the  diode  sampling  the  direct  laser  output  so  the  signal  amplitudes  could 
be  quantitatively  compared.  To  compare  the  signal  wave-forms,  significantly  long 
external  cavities  (>  30  cm)  were  chosen  so  that  the  fundamental  cavity  frequency 
would  be  well  below  the  bandwidth  limit  of  the  500  MHz  oscilloscope. 

Figure  24  shows  two  traces  from  the  RF  signal  analyzer  for  the  feedback  from 
a  normal  mirror  and  the  laser  output.  The  decrease  in  magnitude  of  the  peaks 
above  800  MHz  is  due  to  the  decrease  in  sensitivity  of  the  PIN  diodes  used  in  the 
experiment.  It  can  be  clearly  seen  in  this  figure  that  the  frequency  content  of  both 
the  normal  feedback  and  laser  output  signals  are  the  same.  The  external  cavity  mode 
of  230  MHz  {Lext  =  65  cm)  and  higher  harmonics  are  easily  identified  in  both  signals, 
indicating  that  amplitude  modulation  from  the  onset  of  chaotic  activity  in  the  laser 
is  contained  in  both  the  laser  output  and  the  normal  feedback  signal.  For  a  normal 
dielectric  mirror,  this  response  is  expected  since  these  mirrors  can  be  considered  to 
react  instantaneously  and  have  a  flat  broadband  RF  response. 

The  time  delay  between  the  two  PIN  diode  detectors  is  shown  in  Figure  25  as 
a  correlation  between  the  two  signals.  The  correlation  trace  shows  that  the  feedback 
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Figure  24.  RF  spectra  of  the  feedback  signal  from  a  normal  mirror  (above)  and  the 
laser  output  (below).  Both  signals  contain  the  same  frequency  informa¬ 
tion  indicating  that  there  is  no  filtering  being  done  by  the  normal  mirror. 
The  speed  of  a  normal  mirror  can  be  considered  instantaneous. 
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signal  (middle  trace)  is  highly  correlated  with  the  original  laser  output  (top  trace) 
at  a  point  when  the  feedback  waveform  is  delayed  by  11  ns.  The  calculated  path 
difference  (both  optical  and  electrical)  between  the  two  traces  is  10.9  ns  indicating 
that  all  delays  have  been  properly  accounted. 


At=  11  nS 
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Figure  25.  The  correlation  between  the  two  PIN  diodes  sampling  the  laser  output 
and  normal  feedback  signal  shows  a  peak  at  a  delay  of  11  ns.  This  delay 
corresponds  exactly  to  the  optical  and  electrical  path  differences  between 
the  two  detectors. 


Figure  26  shows  a  similar  RF  spectrum  for  the  case  of  the  CAT  PCM  for  an 
external  cavity  of  87  cm.  The  small  overall  magnitude  of  the  two  traces  is  a  result 
of  the  lower  intensity  necessary  by  the  CAT  PCM  to  undamp  the  external  cavity 
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frequency.  Since  the  phase  conjugate  feedback  couples  easier  into  the  laser  cavity 
than  the  normal  mirror,  significantly  lower  feedback  levels  were  used  to  prevent 
the  laser  from  going  fully  chaotic.  Therefore  the  ND  filter  in  front  of  the  laser 
output  detector  was  increased  to  balance  the  two  signals  for  comparison  purposes. 
Regardless  of  the  absolute  magnitude  of  the  two  signals,  it  is  easily  seen  that  the 
PCM  feedback  signal  also  contains  the  fundamental  external  cavity  frequency  plus 
several  harmonics.  Figure  27  shows  the  correlation  (bottom  trace)  between  the  direct 
laser  output  (top  trace)  and  the  CAT  phase  conjugate  signal  (middle  trace).  As  was 
shown  for  the  normal  dielectric  mirror,  the  two  signals  have  a  strong  correlation 
peak  corresponding  to  the  path  differences  between  the  two  PIN  detectors.  In  other 
words,  the  CAT  PCM  feedback  signal  contains  the  RF  frequency  information  and 
correlates  well  with  the  direct  laser  output  and  therefore,  performs  on  the  order  of 
the  normal  dielectric  mirror. 

5.2.3  External  Cavity  Length  and  Frequency  Dependency.  One  of  the 
easiest  parameters  to  change  in  the  set  of  nonlinear  equations  is  the  time  delay 
due  to  the  position  of  the  feedback  mirror  in  the  external  cavity.  According  to  the 
analysis  of  Auyeung  et  al,  the  mode  structure  of  the  external  cavity  bounded  by  a 
PCM  depends  on  whether  the  PCM  is  established  by  degenerate  or  non-degenerate 
four-wave  nonlinear  interaction  (4).  For  the  degenerate  case,  such  as  the  self-pumped 
CAT  and  DPCM  configurations  used  in  this  research,  the  external  cavity  will  support 
a  frequency  spacing  of  c/{2Lext)-  This  is  obtained  by  requiring  a  single-round  trip 
self-consistent  solution  for  the  complex  radius  of  curvature  of  the  wavefront.  This 
spacing  is  also  the  same  for  a  cavity  bounded  by  a  normal  mirror.  In  contrast,  for  a 
non-generate  PCM  the  field  must  take  two  round  trips  in  the  external  cavity  in  order 
to  satisfy  the  self-consistent  solution,  thus  the  mode  spacing  for  the  non-degenerate 
PCM  becomes  c/{ALext)-  This  change  from  the  non-degenerate  case  results  from  a 
frequency  shift  in  the  conjugate  beam  and  passing  through  the  PCM  a  second  time 
cancels  this  effect.  The  non-degenerate  PCM  was  not  investigated  in  this  research. 
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1  div=  180  MHz 

Figure  26.  RF  spectra  of  the  feedback  signal  from  the  CAT  PCM  mirror  and  the 
laser  output.  Both  signals  contain  the  same  frequency  information  indi¬ 
cating  that  there  is  no  filtering  being  done  by  the  CAT  mirror. 
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Figure  27.  Correlation  function  of  the  CAT  PCM  signal  with  the  direct  laser  output 

The  graphs  in  Figures  28  and  29  show  the  effect  of  the  external  cavity  frequency, 
by  varying  the  external  cavity  length  using  the  CAT  phase  conjugate  mirror 
and  a  normal  mirror.  Each  graph  shows  the  theoretical  limit  of  Vext  =  c/ {2Lext)  for 
the  degenerate  PCM  and  normal  mirror,  where  Lext  is  the  optical  path  measured 
from  the  edge  of  the  laser  facet  to  the  entrance  face  of  the  crystal. 

5.2.4  External  Cavity  Resonance  and  Frequency  Locking.  Based  on  the 
mathematical  analysis  presented  in  Section  3.3,  resonance  enhancement  in  the  feed¬ 
back  levels  for  the  first  Hopf  bifurcation  was  expected  as  the  external  cavity  length 
was  extended.  The  external  cavity  was  scanned  using  a  micrometer  translation 
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Ext  Cavity  Freq  (MHz) 


Figure  28.  Plot  of  experimentally  measured  external  cavity  frequency  using  a  de¬ 
generate  CAT  PCM.  The  solid  line  represents  the  theoretical  value  of 
c/ {2Lext)  for  a  cavity  bounded  by  a  degenerate  phase  conjugate  mirror. 
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Figure  29.  Plot  of  experimentally  measured  external  cavity  frequencies  versus  cavity 
length  using  a  normal  mirror.  The  solid  line  represents  the  theoretical 
value  of  c/{2\jext)  for  an  external  cavity  with  a  normal  mirror. 
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stage  and  a  Klinger  translation  stage.  The  micrometer  stage  was  restricted  by  a  2 
cm  throw  which  caused  a  slight  error  in  the  measurement  by  resetting  the  stage  for 
subsequent  scans.  This  problem  was  corrected  by  replacing  the  micrometer  with  a 
Klinger  translation  stage  which  had  a  significantly  larger  throw.  Figure  30  shows 
the  experimental  equipment  used  to  investigate  the  cavity  resonance  using  the  CAT 
geometry.  The  cavity  was  extended  over  8  cm,  in  1  mm  increments,  to  capture  at 
least  two  resonance  cycles  predicted  by  Eq.  (66).  The  Hopf  point  was  determined  by 
monitoring  the  F-P  display  and  recording  the  feedback  power  when  the  relaxation 
frequency  appeared.  The  phase  conjugate  feedback  was  then  allowed  to  grow  until 
the  external  cavity  frequency  could  be  recorded. 

Figure  31  shows  the  results  of  the  data  collection.  The  top  curve,  which  cor¬ 
responds  to  the  right  vertical  axis,  plots  the  ratio  of  the  relaxation  frequency  to  the 
external  cavity  frequency,  vp_/vext-  The  stair  stepping  nature  of  the  curve  indicates 
frequency  locking  of  ur,  with  the  plateaus  corresponding  to  integer  values  6,  7,  and 
8.  The  lower  part  of  the  curve  shows  the  power  of  the  feedback  signal  at  the  first 
Hopf  bifurcation  when  vr  becomes  undamped.  Resonances  can  be  seen  for  cavity 
lengths  of  22.8  and  26.2  cm  with  the  width  of  the  locking  region  being  3.4  cm.  The 
peak  shown  at  20.2  cm  is  the  trailing  edge  of  an  earlier  peak. 

The  mechanism  allowing  the  integer  steps  in  the  ratio  of  the  two  frequencies 
is  the  change  in  ur.  The  lower  half  of  Figure  32  presents  the  variation  in  ur  and 
compares  it  to  the  locking  ratio  presented  in  the  previous  figure.  The  solitary  laser 
relaxation  oscillation  frequency  of  3.7  GHz  falls  midway  between  the  maximum  and 
minimum  values.  At  the  points  of  transition  to  a  higher  integer  ratio,  pr  will  attain 
its  maximum  value  and  Apr  =  500  MHz  over  the  full  locking  range. 

5.2.5  Verification  of  the  Conjugate  Signal.  As  mentioned  in  Section  5.2.3, 
the  mode  spacing  for  the  degenerate  PCM  and  the  normal  cavity  is  the  same  c/ {2Lext) 
relationship.  Since  the  mode  spacing  is  the  same,  verification  of  a  phase  conjugate 
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Figure  30.  Experimental  setup  used  to  capture  the  feedback  resonance  peaks  and 
frequency  locking  from  the  CAT  PCM.  The  Klinger  translation  stage  was 
used  to  extend  the  cavity  in  a  span  of  8  cm  in  1  mm  increments. 


73 


Conjugate  Pwr  (uW) 


Relaxation  Freq  (GHz) 


External  Cavity  Length  (cm) 


Figure  32.  Experimental  data  showing  the  frequency  locking  (top  curve,  right  axis) 
of  the  external  cavity  and  relaxation  frequencies  from  the  CAT  PCM. 
The  bottom  curve  shows  the  changes  in  the  relaxation  oscillation  as  the 
external  cavity  is  varied.  The  solitary  relaxation  frequency  of  3.7  GHz 
falls  halfway  between  the  maximum  and  minimum  values. 
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signal  was  necessary  to  rule  out  the  possibility  that  the  feedback  beam  for  the  CAT 
mirror  might  contain  a  normal  reflection  component,  possibly  from  total  internal 
reflection  inside  the  crystal.  Therefore,  several  experiments  were  conducted  to  test 
the  validity  of  the  phase  conjugate  signal  using  the  commonly  accepted  method  of 
image  reconstruction  through  a  phase  distorting  medium  in  the  optical  path  (6). 

To  capture  the  image  of  the  conjugate  beam,  a  BeamCode  camera  replaced  the 
Newport  meter  sampling  the  conjugate  beam.  The  camera  was  also  placed  in  the 
beam  entering  the  crystal  to  record  the  profile  of  the  original  beam.  Figure  33  shows 
the  contour  profile  of  the  original  beam  entering  the  crystal.  The  vertical  fringes 


Figure  33.  Profile  of  the  beam  entering  the  BaTiOs  crystal. 


shown  in  the  figure  are  a  result  of  interference  effects  from  several  ND  filters  placed 
in  front  of  the  camera.  The  basic  profile  of  the  unaltered  beam  is  a  2-D  gaussian, 
with  the  horizontal  direction  having  the  smaller  waist. 

Figure  34  shows  the  conjugate  return  profile  for  several  levels  of  feedback.  In 
(a),  the  relaxation  frequency  has  begun  to  undamp  and  the  conjugate  begins  to 
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Figure  34.  Phase  conjugate  intensity  profiles  captured  using  BeamCode  as  chaos 
develops,  (a)  Relaxation  oscillations  are  undamped,  (b)  external  cavity 
frequency  undamped,  (c)  several  external  cavity  frequencies  undamped, 
(d)  full  chaos  behavior. 
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build  around  the  central  portion  of  the  beam  profile.  The  reconstruction  of  just  the 
central  portion  of  the  original  profile  is  a  consequence  of  the  low  feedback  power 
needed  to  undamp  the  relaxation  oscillation.  With  the  self-pumped  geometry,  low 
feedback  power  is  generated  by  low  input  intensity  to  the  BaTiOs  crystal,  hence 
only  the  most  intense  portions  of  the  image  profile  (central  portion)  are  efficiently 
conjugated.  Increasing  the  feedback  until  the  external  cavity  mode  appears,  the 
conjugate  profile  shown  in  (b)  continues  to  grow  around  the  central  portion  of  the 
profile.  Image  (c)  shows  the  beam  profile  as  the  feedback  is  increased  to  the  level 
where  several  harmonics  of  the  external  cavity  frequency  appear.  When  full  chaos 
in  the  laser  output  is  reached,  (d),  the  conjugate  signal  is  well  defined  and  resembles 
the  original  beam.  This  profile  clearly  shows  the  full  reconstruction  of  the  original 
image. 

The  same  sequence  of  events  was  repeated  with  a  normal  cavity  for  comparison. 
Figure  35  shows  the  beam  profiles  corresponding  to  the  events  of  the  relaxation 
oscillations  undamped  (a),  external  cavity  frequency  undamped  (b),  and  full  chaos 
in  the  laser  output  (c).  As  expected,  the  full  beam  profile  is  clearly  visible  in  all 
cases  since  the  normal  mirror  has  a  flat  intensity  response  to  its  input  beam.  The 
only  difference  between  the  images  in  the  figure  is  the  overall  intensity  of  the  return 
beam  entering  the  camera,  which  increases  from  (a)  to  (c). 

Using  a  piece  of  plastic  as  a  phase  distorting  medium,  the  two  tests  were 
repeated.  Figure  36  shows  a  comparison  of  the  conjugate  and  normal  feedback 
profiles.  In  both  cases  the  images  shown  were  taken  when  the  relaxation  and  external 
cavity  frequencies  have  undamped.  It  was  also  noted  that  the  phase  conjugate  mirror 
with  the  distorter  required  about  the  same  amount  of  feedback  power  to  achieve  the 
undamping  of  the  two  frequencies  as  compared  to  the  PCM  without  the  distorter. 
In  contrast,  the  normal  mirror  required  over  50%  more  feedback  power  to  achieve 
the  same  result. 
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(a)  (b) 


(c) 

Figure  35.  Normal  mirror  profiles  captured  using  BeamCode  as  chaos  develops,  (a) 
Relaxation  oscillations  are  undamped,  (b)  external  cavity  frequency  has 
undamped,  (c)  full  chaos.. 
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(a)  (b) 


Figure  36.  Phase  conjugate  signal  (a)  and  normal  signal  (b)  with  a  phase  distorting 
medium  placed  in  the  beam  entering  the  BaTiOa  crystal.  In  both  cases 
the  relaxation  and  external  cavity  frequencies  have  undamped. 
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The  changes  in  feedback  power  for  the  normal  mirror  is  easily  attributed  to 
scattering  of  photons  off  the  distorting  medium.  Coupling  into  the  laser  is  reduced 
since  the  phase  front  is  badly  distorted,  therefore,  more  feedback  power  is  required 
to  overcome  the  losses  in  the  external  cavity.  The  behavior  of  the  conjugate  signal 
demonstrates  that  the  PCM  is  producing  a  conjugate  beam  which  feeds  back  the 
same  amount  of  power  even  with  the  distorting  medium  in  the  beam  path. 

5.3  Double  Phase  Conjugate  Mirror  Measurements 

Most  of  the  tests  described  in  Section  5.2  for  the  CAT  geometry  were  dupli¬ 
cated  with  the  double  phase  conjugate  mirror  (DPCM).  Figure  37  shows  the  basic 
experimental  setup  used  to  conduct  the  tests.  The  mirror  was  constructed  using  the 
transmission  grating  DPCM  with  input  angles  of  120°  degrees,  with  respect  to  the 
-hc-axis.  A  variety  of  input  angles  were  tried,  most  arrangements  produced  little 
or  no  feedback.  It  appeared  that  120°  was  the  most  optimal  choice  of  angle.  Even 
at  this  setting,  the  maximum  amount  of  feedback  was  small  compared  to  the  CAT 
geometry.  As  a  result,  full  chaotic  behavior  was  never  achieved  with  the  DPCM, 
although  sufficient  feedback  was  obtained  to  operate  in  Levels  I,  II  and  III.  Due  to 
complexities  of  the  DPCM  geometry,  the  external  cavity  was  not  varied  as  exten¬ 
sively  as  the  CAT  geometry,  but  was  varied  enough  to  verify  the  trends  shown  in 
the  CAT  external  cavity  experiment. 

The  experimental  results  of  the  conjugate  feedback  of  the  DPCM  configuration 
show  similar  behavior  to  the  CAT  PCM.  Figure  38  shows  the  external  cavity  mode 
spacing  as  a  function  of  the  round  trip  external  cavity  length.  Measurement  of  the 
cavity  length  is  slightly  different  than  the  CAT  PCM  geometry  since  the  DPCM  is 
an  asymmetric  cavity.  For  the  DPCM,  photons  from  one  beam  are  scattered  off  the 
holographic  grating  in  the  direction  of  the  other  input  beam’s  conjugate  signal,  as 
was  shown  in  Figure  7.  Therefore,  the  external  cavity  is  in  a  ring  formation  as  the 
photons  traverse  through  the  BaTiOa  crystal,  retracing  the  opposite  input  beam’s 
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Figure  37.  Experimental  configuration  used  in  the  double  phase  conjugate  mirror 
experiments. 
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phase  front.  This  is  different  from  the  CAT  geometry  where  the  input  and  conjugate 
beams  share  the  same  path.  Hence  the  round  trip  external  cavity  length  of  2Lext  as 
measured  in  the  CAT  PCM  experiments  now  becomes  L^ing,  where  Lring  is  the  total 
optical  path  length  of  the  beam  leaving  the  laser,  passing  through  the  crystal  and 
then  re-entering  the  laser.  So  for  the  degenerate  DPCM,  the  theoretical  external 
cavity  frequency  spacing  is  restated  as  c/Lring,  as  shown  in  Figure  38. 


Figure  38.  Plot  of  experimentally  measured  external  cavity  frequencies  versus  cavity 
length  using  a  DPCM  mirror.  The  solid  line  represents  the  theoretical 
mode  spacing  value  of  c/Lring  for  a  degenerate  PCM. 


The  RF  response  of  the  DPCM  was  also  tested  in  a  similar  manner  to  the  CAT 
and  normal  mirror  tests  presented  earlier.  Figure  39  shows  the  two  RF  spectra  of 
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Figure  39.  RF  spectra  of  the  feedback  signal  from  the  DPCM  mirror  and  the  laser 
output.  Both  signals  contain  the  same  frequency  information  indicating 
that  there  is  no  filtering  being  done  by  the  DPCM  mirror 

the  DPCM  conjugate  beam  and  the  direct  laser  output.  Since  the  conjugate  levels 
were  extremely  low,  the  ND  filter  in  front  of  the  PIN  diode  sampling  the  laser  output 
was  removed  to  identify  the  RF  signals  clearly.  Consequently,  no  comparison  can  be 
made  between  the  amplitudes  of  the  two  spectra  shown  in  the  figure.  However,  the 
spectra  do  show  the  primary  peak  of  Uext  =  350  MHz  which  corresponds  to  Lring  = 
85.7  cm.  A  second  weak  peak  in  the  DPCM  spectrum  can  also  be  seen  at  2Pext  =  700 
MHz.  No  other  peaks  can  be  seen,  but  the  sensitivity  of  the  detector  dramatically 
decreases  for  frequencies  >  800  MHz,  which  would  explain  the  lack  of  peaks  beyond 
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2ug^t.  The  direct  laser  output  spectrum  shows  this  dramatic  decrease  in  sensitivity 
for  u  =  1400  MHz. 

The  next  chapter  takes  the  raw  experimental  data  and  makes  a  comparison 
against  the  predictions  of  the  analytical  solutions  and  computer  models  discussed  in 
Chapters  III  and  IV. 
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VI.  Experimental  and  Modeling  Results 

The  culmination  of  this  research  was  to  provide  insight  into  the  processes  of 
chaotic  development  in  laser  diodes  with  phase  conjugate  feedback  and  the  validation 
of  a  computer  model  that  could  accurately  predict  this  nonlinear  behavior.  This 
chapter  explains  the  improvement  made  to  the  LKCON2.C  model  and  compares  the 
experimental  data  with  the  model  outputs. 

6.1  Gain  Saturation  Improvement 

As  mentioned  in  Section  4.2,  an  inaccuracy  was  identified  with  the  original 
computer  model  for  long  external  cavities.  In  this  region,  the  model’s  prediction 
of  the  external  cavity  frequency  quickly  diverged  from  experimental  data.  The  re¬ 
sult  for  long  cavities  {L^xt  ~  15-30  cm)  was  an  external  cavity  frequency  spacing 
that  was  underestimated  by  as  much  as  one-fourth  of  the  measured  experimental 
frequency  spacing.  The  same  deficiency  may  also  exist  for  short  cavities;  however, 
at  short  cavity  lengths,  the  dominant  route  to  chaos  is  through  period  doubling  of 
the  relaxation  oscillation  which  is  seen  in  the  model.  Simulations  with  short  cavi¬ 
ties  show  that  the  period  doubling  of  the  relaxation  oscillation  may  actually  mask 
the  incorrect  external  cavity  frequency  since  Level  IV  behavior  is  reached  before  the 
external  cavity  frequency  undamps.  Since  no  experimental  data  had  been  collected 
for  Lext  <  15  cm,  validation  of  the  external  cavity  frequency  spacing  below  15  cm 
was  not  possible. 

The  original  computer  model  was  based  on  an  underlying  assumption  in  the  set 
of  coupled  nonlinear  equations  that  the  gain  coefficient  presented  in  Eq.  (18)  varied 
linearly  with  AV.  This  first  order  linear  approximation  was  an  insufficient  estimate 
of  the  gain  coefficient.  A  more  realistic  representation  requires  the  inclusion  of  a 
small  signal  gain  saturation  term  Y  to  limit  the  effect  of  the  electron  hole  number 
on  the  photon  number  as  shown  below. 
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The  inclusion  of  the  gain  saturation  term  only  affects  da/ds  and  dN/ds  pre¬ 
sented  in  Eqs.  (52)-(54).  The  corrected  set  of  nonlinear  equations  used  in  the 
LKCON2.C  model  are 


da  1  2N  +  1  ^ 

ds  2  yi  +  YV{a^  -1)  \ 

-I-  ria{s  -  r)  cos{(/)(s  -  r)  -I-  (j){s)  -  ^pcm  -  2Alj{s  -  r/2)) 


(71) 


and 


T"  =  r. 

ds 


N- 


(1  +  2iV) 
yi  +  YV(a'‘  -  1) 


(72) 


where  Y  is  the  dimensionless  small  signal  gain  saturation  term  {Y  >0).  Setting 
Y  =  0.0  recovers  the  original  set  of  equations.  Notice  that  for  F  >  0,  the  amplitude  of 
the  field  reduces  the  magnitude  of  the  stimulated  gain  term  included  in  N.  Increasing 
the  value  of  Y  causes  the  stimulated  gain  to  saturate  more  quickly. 


Simpson  et  al  measured  the  small  signal  gain  term  to  be  F  =  0.017  in  a  solitary 
laser  diode  similar  to  the  SDL  5412  laser  diode  used  in  our  experiments  (24).  It  was 
found  through  comparison  of  the  model  and  experimental  data  that  a  value  of  F  = 
0.017  was  adequate  for  the  SDL  5412.  Even  though  it  was  beyond  the  scope  of  this 
research  to  specifically  measure  the  small  signal  gain  term,  the  excellent  agreement 
between  the  computer  model  results  and  experimental  data  confirms  that  the  value 
of  F  =  0.017  is  reasonable. 


With  the  gain  saturation  term  added  to  the  model  a  comparison  between  ex¬ 
perimental  data  and  model  predictions  can  be  accomplished.  The  next  several  sec¬ 
tions  present  the  analysis  of  the  model  and  experimental  data  for  the  external  cavity 
frequency  spacing,  changes  in  relaxation  oscillation,  locking  ratio  and  feedback  res¬ 
onance.  A  qualitative  comparison  between  the  bifurcation  map  and  experimental 
data  is  discussed  in  the  last  section. 
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6.2  External  Cavity 

Figure  40  shows  the  experimental  (closed  diamonds)  and  model  data  (dotted 
line)  for  the  external  cavity  frequency  as  a  function  of  the  external  cavity  length  in 
the  range  Lext  =  20-27  cm.  Also  shown  in  the  figure  is  the  v  =  cj {2Lext)  hue  which 


Y=  0.017 

Figure  40.  Comparison  between  the  model  (open  diamond)  and  experimental  (dot¬ 
ted  line)  data  for  the  external  cavity  frequency,  v^xt^  as  a  function  of  the 
external  cavity  length,  Lext- 


is  the  theoretical  mode  spacing  for  degenerate  phase  conjugate  feedback  discussed  by 
Auyeung  et  al.  As  shown  in  the  figure,  there  is  excellent  overall  agreement  between 
the  experimental  and  model  results.  The  maximum  difference  between  these  two 
data  sets  is  on  the  order  of  25  MHz  which  is  outside  the  3  MHz  resolution  of  the 
spectrum  analyzer.  The  additional  error  between  the  data  sets  is  attributed  to  the 
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limited  accuracy  of  the  FFT  in  the  model  (±  20  MHz)  and  to  a  lesser  extent,  to  an 
error  in  measuring  the  experimental  external  cavity  (±  2  mm). 

Despite  the  small  differences  between  the  experimental  and  model  data  sets, 
there  still  remains  an  overall  downward  shift  of  the  frequencies  relative  to  the  c/ {2Lext) 
line.  The  downward  shift  in  the  frequency  data  cannot  be  attributed  to  a  frequency 
bias  of  the  RF  spectrum  analyzer.  The  analyzer  was  tested  by  splitting  the  soli¬ 
tary  laser  output  and  up-shifting  one  of  the  beam’s  frequency  by  48  MHz  using  an 
acousto-optic  modulator.  The  two  beams  were  then  recombined  onto  the  PIN  diode 
detector  and  the  beat  signal  between  the  two  laser  beams  was  displayed  on  the  RF 
spectrum  analyzer.  The  manufacturer’s  specification  for  the  instrument  is  a  display 
resolution  <  3  MHz  over  the  full  span  (DC-1.8  GHz)  of  the  analyzer.  The  48  MHz 
beat  signal  was  detected  well  within  this  specification. 

In  addition,  a  horizontal  translation  of  the  experimental  and  model  data  to 
coincide  with  the  c/{2Lext)  line  is  not  reasonable.  Although  a  cavity  length  mea¬ 
surement  error  of  ±  2  mm  does  exist,  this  error  is  well  below  the  >3  cm  translation 
needed  to  move  the  experimental  data  onto  the  c/{2Lext)  line.  The  most  reasonable 
explanation  for  the  overall  downward  shift  of  both  the  experimental  and  computer 
model  data  is  the  presence  of  other  non-linear  processes,  such  as  gain  saturation, 
that  have  not  been  included  in  Auyeung’s  analysis.  The  development  of  the  c/ {2Lext) 
line  was  based  on  the  self-consistent  solution  of  the  radius  of  curvature  of  the  laser 
beam  in  a  simple  cavity  formed  between  the  phase  conjugate  mirror  and  the  gain 
medium.  No  other  processes  or  limitations  on  the  lasing  medium  or  the  phase  conju¬ 
gate  mirror  were  considered.  Furthermore,  this  downward  shift  in  was  also  seen 
in  the  experimental  results  from  the  DPCM  and  to  a  lesser  degree  in  the  normal 
mirror  experiments.  This  suggests  the  shift  is  strictly  a  laser  and  external  cavity 
phenomenon  and  is  independent  of  the  feedback  method. 
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6.3  Relaxation  Oscillation  and  Locking  Ratio 

Figure  41  shows  the  comparison  of  the  raw  experimental  data  and  the  LK- 
CON2.C  model  for  changes  in  the  relaxation  oscillation  as  the  external  cavity  length 
is  varied.  The  relevant  model  parameters  are  listed  at  the  bottom  of  the  figure.  Both 
curves  have  the  same  ALext  —  3-5  cm  between  successive  transitions,  as  well  as  Aur 
=  500  MHz  as  the  max-min  transition  of  the  relaxation  oscillation.  The  general 
features  of  the  experimental  data  seem  to  be  captured  by  the  model  curve.  However, 
it  is  apparent  from  the  figure  that  there  seems  to  be  a  shift  of  the  model  data  both 
in  frequency  and  external  cavity  length  relative  to  the  experimental  data. 

Earlier  in  Section  3.3  it  was  shown  that  the  Hopf  frequency,  CIr,  emerging  out 
of  the  first  Hopf  bifurcation  varied  around  the  solitary  relaxation  value.  As  shown 
in  Figure  41,  this  corresponds  to  a  =  4.28  GHz  for  the  model  and  3.72  GHz  for 
the  experimental  data.  As  a  cross-check,  the  solitary  relaxation  frequency  can  also 
be  extracted  from  the  span  between  adjacent  transitions.  Since  the  span  between 
resonant  transitions  is  27r,  an  expected  value  for  the  experimental  ur  using  Eq.  (68) 
is  calculated  to  be: 

i^R  =  =  4-28  GHz  (73) 

lIALext 

This  is  obviously  different  than  the  measured  value  of  vr  =  3.72  GHz.  Eigure  42 
shows  a  diagram  of  the  E-P  display  for  an  arbitrary  cavity  length.  Erom  the  display, 
it  is  not  possible  to  determine  which  frequency  is  to  be  associated  with  adjacent  E-P 
cavity  modes.  Since  the  FSR  of  the  F-P  is  8  GHz,  the  natural  association  is  to  assume 
that  the  two  frequencies  displayed  are  at  3.72  GHz.  However,  the  model  provides 
insight  into  the  proper  interpretation  of  the  F-P  display,  which  is  overlapping  modes 
with  frequencies  at  4.28  GHz.  Therefore  the  data  must  be  reinterpreted  with  a 
solitary  relaxation  frequency  of  4.28  GHz,  which  also  agrees  with  the  calculated 
value  based  on  the  span  between  resonant  transitions. 
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3.72  GHz 


4.28  GHz 


Figure  42.  Diagram  of  a  F-P  display  showing  an  8  GHz  FSR  between  adjacent  F-P 
modes  and  two  frequencies.  From  the  F-P  display  it  is  not  possible  to 
determine  if  the  two  frequencies  are  at  3.72  GHz  or  4.28  GHz. 
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The  reinterpretation  of  the  experimental  data  accounts  for  the  vertical  shift 
between  the  computer  model  and  the  experimental  data.  To  address  the  horizontal 
differences,  we  need  to  apply  the  approximation  of  a  — >  oo  made  in  Section  3.3. 
It  was  shown  that  this  approximation  simplifies  the  equations  for  the  relaxation 
frequency  emerging  out  of  the  first  Hopf  bifurcation.  Accordingly,  the  a  parameter 
in  the  model  was  increased  to  a  relatively  high  value  of  10.  The  effect  of  the  increase 
in  a  was  to  shift  the  model  data  horizontally  along  Lext  so  that  the  transition  points 
at  22.7  and  26.2  cm  match.  In  effect,  the  increase  in  a  reduces  the  dependency  of 
Eq.  (67)  on  the  cos('d)  term  in  the  denominator  and  hence  transitions  in  Q.r  (or  vr) 
will  occur  when  VLt  =  {2n  +  l)7r.  Increasing  a  beyond  a  value  of  10  in  the  model  has 
no  significant  contribution  to  the  shift  in  the  data  since  the  denominator  of  Eq.  (67) 
goes  as  arctan(a).  Figure  43  shows  the  result  of  incorporating  both  the  increase 
of  a  in  the  model  data  and  the  reinterpreted  experimental  values.  The  excellent 
agreement  between  the  two  curves  is  evident  in  the  figure. 

Figure  44  shows  the  integer  values  of  the  locking  ratio  for  the  experimental  and 
model  data.  The  locking  ratio  values  for  the  model  data  have  been  calculated  using 
the  a  =  10  model  data  and  the  ratios  for  the  experimental  data  have  been  calculated 
using  the  relaxation  frequencies  shown  in  Figure  43.  As  shown  in  Figure  44,  the  graph 
of  the  experimental  data  has  a  whole  integer  trend  with  values  of  7,  8,  and  9  over  the 
range  of  =  20-28  cm  while  the  model  data  has  half-integer  values  of  6.5,  7.5  and 
8.5  over  the  same  range.  Earlier  in  the  external  cavity  frequency  comparison,  the 
experimental  data  had  been  shown  to  be  on  the  order  of  25  MHz  below  the  model 
data.  This  small  difference  is  significant  enough  to  increase  the  locking  ratio  of  the 
experimental  data  by  a  half-integer.  Re-scaling  the  experimental  data  to  coincide 
with  the  model  trend  would  lower  the  experimental  locking  ratio  values  by  a  half 
integer.  Therefore,  the  disparity  in  the  locking  ratios  is  attributed  to  the  small,  but 
consistent,  differences  in  the  external  cavity  frequency  between  the  model  and  the 
experimental  data. 
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»'r(GHz] 


^  External  Cavity  Length  (cm) 

a  =  10.0 

T  =  316.85 
Y=  0.017 


Figure  43.  Comparison  between  experimental  and  model  data  for  the  changes  in 
the  relaxation  oscillation  frequency  using  the  increased  value  of  a  in  the 
model  and  reinterpretation  of  the  experimentally  measured  pr. 
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6-4-  Feedback  Resonance 

The  last  comparison  to  be  discussed  is  between  the  experimental  external  cav¬ 
ity  resonance  data  and  Eq.  (69)  describing  r]H  using  the  approximation  a  — ^  oo. 
Figure  45  compares  the  measured  experimental  and  model  results  for  external  cavity 
lengths  of  20-28  cm.  The  experimental  phase  conjugate  feedback  power  was  nor¬ 
malized  using  Eq.  (48)  for  k.  The  laser  cavity  coupling  efficiency  in  Eq.  (48)  was 
assumed  to  be  =  1  and  a  value  for  the  laser  facet  reflectivity,  R2  =  0.3,  was  chosen 
based  on  manufacturer  specifications.  The  resultant  value  for  k  was  then  rescaled 
using  the  scaling  presented  in  Section  3.1  to  convert  the  experimentally  measured 
feedback  power  into  the  dimensionless  parameter  77. 

The  significant  attribute  of  this  figure  is  the  overall  fit  of  the  periodicity  be¬ 
tween  the  data  sets.  Qualitatively,  Eq.  (69)  accurately  captures  the  general  trend 
of  the  data.  Again,  the  3.5  cm  peak-to-peak  resonance  spacing  in  the  experimental 
data  suggests  a  calculated  vr  —  4.28  GHz.  The  model  representation  shown  in  the 
figure  was  calculated  using  ur  =  4.28  GHz  which  yields  an  excellent  fit  to  the  data. 
This  value  for  the  relaxation  frequency  also  agrees  with  the  value  established  earlier 
in  Section  6.3  when  the  LKCON2.C  model  outputs  for  the  relaxation  frequency  were 
compared  to  the  experimental  values.  For  comparison,  using  a  relaxation  oscillation 
value  Ur  =  3.72  GHz  in  Eq.  (69)  yields  a  resonance  spacing  of  4  cm,  which  was  orig¬ 
inally  shown  in  Figure  10  when  the  equation  for  rjR  was  first  introduced.  A  curve 
fit  to  the  experimental  data  with  a  resonance  spacing  of  4  cm  is  poor  since  only  one 
peak  can  be  matched  with  that  spacing. 

A  notable  conclusion  from  this  graph  is  that  the  experimental  data,  the  LK- 
CON2.C  model,  and  the  approximate  solution  for  rjH  given  in  Eq.  (69)  all  agree  over 
this  range.  Since  Eq.  (69)  offers  a  convenient  and  compact  form  to  investigate  cavity 
resonance  behavior,  a  validating  test  would  be  to  apply  the  external  cavity  plot  to 
a  second  set  of  data  and  look  for  resonance  point  matches.  As  mentioned  in  the 
experimental  descriptions,  external  cavity  resonance  data  were  also  taken  indepen- 
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dently  using  a  manual  translation  stage  which  only  had  a  2  cm  throw.  Figure  46 
plots  the  translation  stage  data  against  Eq.  (69)  for  the  ranges  25-27  cm  and  28-30 
cm  which  accurately  shows  the  correlation  between  the  experimental  data  and  the 
approximate  solution. 

6.5  Bifurcation  Diagram 

The  direct  experimental  construction  of  a  bifurcation  diagram,  depicting  the 
extreme  values  of  the  laser  output  amplitude  versus  feedback  rate,  was  not  possible 
due  to  data  acquisition  equipment  limitations.  However,  a  qualitative  comparison 
between  experimental  data  and  the  LKCONBI2.C  model  was  accomplished  and  nsed 
throughout  the  experimental  work.  The  comparison  uses  experimental  RF  and  F-P 
spectra  to  identify  key  events  in  the  model  bifurcation  map,  such  as  relaxation 
oscillation  undamping,  the  external  cavity  frequency  appearing  and  total  chaos.  The 
experimental  feedback  power  associated  with  these  spectral  events  can  then  be  used 
to  locate  the  position  in  the  model  bifurcation  map.  This  procedure  was  also  used 
in  reverse  where  model  predictions  could  be  used  to  forecast  expected  experimental 
values  at  least  to  within  an  order  of  magnitude. 

Figure  47  shows  the  model  bifurcation  map  for  L^xt  =  20  cm.  The  diamonds 
overlaid  on  top  are  the  locations  of  the  key  events  found  in  the  experimental  spectra. 
The  phase  conjugate  feedback  power  has  been  normalized  in  a  similar  manner  as 
described  in  the  previous  section  so  that  the  feedback  power  is  now  displayed  as  rj. 
Table  8  lists  the  values  of  the  key  events  for  the  experimental  and  model  data.  The 
emergence  of  the  external  cavity  frequency  at  =  0.0065  is  difficult  to  see  in  the 
map  due  to  the  limited  resolution  of  the  graphic,  but  was  confirmed  by  examining 
the  model  RF  spectrum.  No  attempt  was  made  to  calibrate  the  sensitivity  of  the 
F-P  or  spectrum  analyzer  instruments,  therefore,  exact  pinpointing  of  the  bifurcation 
points  in  the  experimental  data  is  not  reasonable.  Nonetheless,  there  is  a  general 
agreement  between  the  experimental  data  and  the  LKCONBI2.C  model. 
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Figure  47.  Bifurcation  map  created  by  the  LKCONBI2.C  model  for  Lext  =  20 
cm  overlaid  with  points  from  the  experimentally  observed  bifurcations 
(closed  diamonds).  The  experimental  bifurcation  points  were  obtained 
by  watching  the  RF  or  F-P  spectra  and  noting  the  conjugate  feedback 
power  as  frequencies  emerged. 
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Table  8.  Identification  of  key  bifurcation  events  sliowii  in  Figure  47. 


Key  Event 

j^model 

^exp 

Relaxation  Oscillation  Emerges 

0.0051 

0.0063 

External  Cavity  Frequency  Emerges 

0.0065 

0.008 

Full  Chaos 

0.013 

0.014 

This  simple  comparison  is  obviously  not  a  complete  validation  of  the  LK- 
CONBI2.C  computer  model.  Nevertheless  it  demonstrates  good  qualitative  agree¬ 
ment  between  predicted  results  and  experimental  data.  In  order  to  evaluate  the 
model  quantitatively,  experimental  reconstruction  of  the  bifurcation  maps  is  essen¬ 
tial. 
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VII.  Conclusions  and  Recommendations 

Presented  in  this  research  is  an  examination  of  the  nonlinear  response  of  a 
laser  diode  subject  to  weak  degenerate  phase  conjugate  feedback.  The  development 
of  chaos  in  the  laser  diode  was  investigated  through  numerical  computer  modeling 
to  predict  the  laser  output  spectra  and  through  experimental  verification  of  the 
computer  models  using  a  BaTiOa  crystal  as  the  phase  conjugate  feedback  mirror. 
Two  methods  of  establishing  the  PCM  were  used  in  the  experiments  to  compare 
the  feedback  effects  due  to  self-pumping  the  BaTiOs  crystal  and  by  double  phase 
conjugation. 

It  was  shown  that  a  phase  conjugate  mirror  in  either  the  self-pumped  or  double 
pumped  configuration  with  BaTiOs  has  a  frequency  response  snfficient  to  reflect  the 
RF  external  cavity  frequencies  contained  in  the  laser  output.  Although  BaTiOs 
is  inherently  a  slow  crystal  in  forming  the  holographic  grating  needed  to  establish 
the  phase  conjugate  mirror,  once  the  grating  is  established  it  can  support  small 
perturbations  in  the  laser  output  as  the  laser  evolves  into  a  chaotic  state.  Hence 
the  BaTiOs  PCM  performs  as  “fast”  as  a  normal  dielectric  mirror.  In  addition,  the 
feedback  parameter  k  that  couples  the  external  cavity  delay  into  the  rate  equation  for 
the  laser  photon  number  can  be  modeled  accurately  as  a  time  independent  parameter. 

The  experimental  data  also  verified  model  results  using  variations  of  the  ex¬ 
ternal  cavity  as  the  control  variable.  Through  the  use  of  experimental  data  it  was 
clear  that  the  model  needed  to  be  modified  to  include  gain  saturation  of  the  laser 
medium  to  show  accurately  the  dependency  of  the  external  cavity  freqnency  on  the 
length  of  the  external  cavity.  Without  including  gain  saturation,  model  predictions 
of  Uext  for  long  external  cavities  grossly  underestimated  values  obtained  experimen¬ 
tally.  Therefore,  the  electron-hole  gain  coefficient,  initially  assumed  to  vary  linearly 
with  AN,  actually  has  a  more  complex  dependency.  Furthermore,  the  inclusion  of 
the  gain  saturation  on  the  electron-hole  population  corrected  the  model  for  other 
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effects  such  as  frequency  locking  between  the  relaxation  oscillation  and  external  cav¬ 
ity  frequency,  changes  in  the  relaxation  oscillation  frequency,  and  feedback  resonance 
of  the  external  cavity.  These  phenomena  were  all  verified  experimentally  using  the 
self-pumped  phase  conjugate  geometry. 

The  experimental  and  analytical  results  presented  in  this  document  provide 
valuable  insight  into  areas  of  future  research.  As  discussed  earlier,  there  is  excel¬ 
lent  agreement  between  model  predictions  and  experimental  data  for  the  case  of 
degenerate  phase  conjugate  feedback.  The  investigation  of  the  onset  of  chaos  from 
non-degenerate  phase  conjugate  feedback  would  be  most  advantageous  for  a  more 
robust  computer  model.  Non-degenerate  phase  conjugate  feedback  can  produce  ef¬ 
fects  on  the  laser  not  possible  with  degenerate  feedback  or  feedback  from  a  normal 
dielectric  mirror.  The  fold  bifurcation,  external  cavity  mode  spacing  of  and 
locking  of  the  laser  frequency  to  the  pump  frequency  are  all  unique  areas  that  need 
experimental  verification  for  the  models. 

Moreover,  it  was  discussed  earlier  that  the  LKCON2.C  model  was  the  only 
model  verified  experimentally  due  to  the  limitation  of  the  data  acquisition  capacity 
in  the  experimental  equipment.  The  LKCONBI2.C  bifurcation  model  should  also 
be  experimentally  verified  for  completeness.  This  will  require  a  data  acquisition  sys¬ 
tem  to  include  an  ultra-fast  photodiode  with  a  cut-off  frequency  >  8  GHz  and  a 
matching  analog-to-digital  (A/D)  digitizer  to  record  the  amplitude  time  series  for 
the  bifurcation  map.  Also,  an  RF  spectrum  analyzer  with  at  least  a  6  GHz  band¬ 
width  is  recommended  to  capture  the  external  cavity  frequency  and  the  relaxation 
oscillation  frequency  in  the  same  spectrum.  This  would  allow  comparison  of  the 
amplitudes  of  the  two  frequencies  in  addition  to  pinpointing  the  level  of  feedback 
when  each  frequency  emerges  while  minimizing  relative  scale  errors  associated  with 
separate  equipment.  These  bifurcation  points  could  not  be  adequately  identified  in 
the  current  experimental  equipment  since  only  i>ext  was  collected  in  the  RF  spec¬ 
tra.  Although  there  was  excellent  agreement  in  the  frequency  content  between  the 
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RF  and  F-P  spectra,  the  sensitivity  and  resolution  of  the  F-P  are  less  than  the  RF 
spectrum  analyzer. 

Two  additional  areas  for  proposed  research  include  the  experimental  investi¬ 
gation  of  short  cavities  {uext  >  ^r)  and  the  experimental  determination  of  the  small 
signal  gain  parameter.  The  investigation  of  short  cavities  poses  a  challenge  to  en¬ 
gineer  optics  into  an  external  cavity  with  lengths  shorter  than  5  cm.  These  short 
external  cavities  can  also  further  validate  the  fast  response  of  the  PCM  by  driving  the 
external  cavity  frequency  higher  than  the  relaxation  oscillation.  If  indeed  the  phase 
conjugate  return  beam  contains  the  same  frequency  content  as  the  laser  output,  even 
for  short  cavities,  then  the  BaTiOs  PCM  behaves  nearly  instantaneously. 

One  of  the  conclusions  drawn  from  this  research  is  that  the  small  gain  satu¬ 
ration  term  plays  an  important  role  in  modeling  nonlinear  behavior  due  to  phase 
conjugate  feedback.  The  value  of  the  gain  saturation  term  used  in  this  research  was 
based  on  previous  work  by  Simpson  et  al  on  a  similar  laser  diode.  A  recommended 
area  of  future  research  is  to  investigate  the  gain  saturation  process  in  a  single  mode 
laser  diode  to  fully  understand  the  nonlinear  effects  on  the  gain  due  to  feedback  and 
improve  the  computer  models. 

In  conclusion,  the  data  presented  in  this  research  has  increased  the  under¬ 
standing  of  the  onset  of  deterministic  chaos  in  a  semiconductor  laser  diode  with 
phase  conjugate  feedback.  As  the  theory  and  computer  models  describing  nonlin¬ 
ear  behavior  matures,  eventual  control  of  the  chaotic  laser  will  be  practical.  Being 
able  to  control  the  output  of  a  chaotic  laser  opens  real  possibilities  for  encoding 
and  encrypting  information  in  the  chaos  which  could  transition  a  once  undesirable 
condition  into  a  useful  product. 
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Appendix  A.  MathCad  Templates 


MathCad  template  to  examine  the  frequency  and  feedback  rate  at  the  Hopf  bifurcation  point 
"Cap  Omega"  is  the  dimensionless  frequency  emerging  from  the  Bifurcation 
and  "omega"  is  the  scaled  frequency  in  GHz.  "Eta  is  the  dimensionless  feedback  rate  at  the 
Hopf  point. 

Constant  Parameters: 

Pumping:  P  =2.0  Speed  of  Light  (cm/sec):  c  =310'° 

Alpha,  a  4.0  Photon  Lifetime  (sec):  ip  -4.16610  '^ 

Carrier  Lifetime  (sec):  x  tsj  =  1 .32- 1  O'  ^ 


Step  Size  of  L  increment:  AL  =  .01 
Max  L:  Lj^^^  31  Min  L:  L^j„  =20 

^ ^  max  ^  min  -i 

Number  of  steps  taken:  N  - - - N  =  1.1’10 


External  Cavity  Length  (cm):  L  =  L  mjn.L  min  +  ^  max 

Array  Index  for  saving:  i  =0..N  1  j  =0..  1 
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L 


For  the  condition  that  alpha  tends  to  infinity: 

/cort(L)\ 


CO  h(L)  :  = 


1  -(-2-P 


•CO  R'tanl 


Feedback  rate  at  the  Hopf  Point: 


i1h(L)  :  =  - 


e(l-t-2P) 


fl  -t-  a^-cos(0  -I-  atan(a))-cos 


U  r-t(L)\ 


cos(  0  )  I 
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Constants:  c  3-10'°  ip  MO 

Laser  Parameters:  Measured  Relaxation  Freq:  v  ^  ^  3.72- 10^ 

Scaled  Relaxation  Freq:  co  ^  =2,-:rt  v p  (0^  =  0.02337 

Alpha:  a  4 
Pumping  Term:  P  2 
External  Cavity  Length  (cm):  L  3.0 

2‘L  gjjj 

Scaled  External  Cavity  Length:  t  ■  -  x  =200 

CT  p 

C  9 

Calculated  External  Cavity  Frequency:  -ji —  ''ext "2.5*  10 

ext 


Range  of  values  for  omega:  (o  =  .02,  .019...02 

|(0| 

Fold  bifurcation  Lines  ri(co)  ~ 

which  show  the  amount 
of  feedback,  eta: 
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Graphical  Solution  to  frequency  and  phase  lock  equation: 
Range  of  values  for  Theta:  0  :=0,.l ..  S-ti 


Let  eta  =  1  il  :=  1 

Phase  lock  equation:  co(6)  -(-a^-sin(6-)-atan(a)) 
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Appendix  B.  LKC0NBI2.  C  Source  Code 

/* - lkconbi2 .  c - 

/* - Performs  bifurcations  for  the - Lang  Kobayashi - 

/* — model  of  semiconductor  equations  with - 

/* - conjugate  feedback - 

/*  and  gain  saturation 

/* — This  is  a  code  for  delay-differential  equations - 

/* - with  the  appropriate  scalings  perfomed - 

/* - 10  Aug  95 - 

/* -  Original  LICONBI.C  program  developed  by  */ 

/* - Dr  Gavrielides  Phillips  Laboratory 

/*  Modified  by  G.  Hengst  5  Dec  96  to  read  in  */ 

/*  default  paramieters  so  that  the  codes  does  not  */ 

/*  have  to  be  re-compiled  on  every  constant  change  */ 

/*  Modified  by  G.  Hengst  22  Apr  97  */ 

/*  to  include  gain  saturation  effects  */ 


- */ 

-*/ 

*/ 

-*/ 

*/ 

-*/ 

-*/ 

-*/ 


#include  <math.h> 
tinclude  <stdio.h> 
#include  <malloc.h> 


/* 

******  this  section  is  now 

'  read  in  from  LKCONBIF.CFG  file 

/* 

#def ine 

BEGIN  0.0  */ 

/* 

#def ine 

END  500.0  */ 

/*  500.0  */ 

/* 

#def ine 

EXTIME  400.0  */ 

/*  400.0  */ 

/* 

#def ine 

NPOINTS  16000  */ 

/*  16000  */ 

/* 

#def ine 

KMAX  1  */ 

/* 

#def ine 

HTAMIN  0.00018 

*/ 

/*  .00003  .0017  */ 

/* 

#def ine 

HTAMAX  .001  */ 

/*  .001  .0006  */ 

/* 

#def ine 

HTASTEP  100  */ 

/* 

#def ine 

tao  1428.57  */ 

/*  4000.0  */ 

/* 

#def ine 

pi  2.33  */ 

/*  .75  */ 

/* 

#def ine 

alfa  4.0  */ 

/*  10.0  */ 

/* 

#def ine 

taoex  1095.23  */ 

/*  49.67  132.460/2.0  */ 

/* 

#def ine 

delta  0.0  */ 

/* 

#def ine 

phipcm  0.0  */ 

/* 

#def ine 

model  1  */ 

/* 

model=0  phase  equation  */ 

*/ 


/*  model=l  full  model 


*/ 


no 


float  BEGIN,  END,  Y,  EXTIME; 
int  NPOINTS,  KMAX; 
float  HTAMIN,  HTAMAX,  Y; 
int  HTASTEP; 

float  tao,  pi,  alfa,  taoex,  delta,  phipcm,  ext_length; 
float  carrier_lif e ,  photon_life; 
int  model,  usernd; 


void 

rkdumb ( ) ; 

void 

rk4(); 

float 

^vector  0 ; 

void 

free_vector() ; 

void 

nrerror 0 ; 

float 

**matrix() ; 

void 

f ree_matrix() ; 

int 

init_stack() ; 

int 

write_stack() ; 

float 

ranlO  ; 

float  **stack, *save ; 
float  lastde=0 . 0 ,taupi ; 
FILE  *outfile; 

FILE  *cfgfile; 

int  nstack, iseed=-872581 ; 

struct  given  { 

int  dim; 
float  out; 

}  data; 


/*  main(argc , argv)  */ 
main(void) 

/*int  argc; 
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char  *argv  []  ;  */ 

{ 

void  derivsO; 
void  exitO  ; 
float  *ystar; 
register  int  i; 
float  hjtaUjtpi; 
int  dimen, j ; 


/* - */ 

/*  This  section  has  been  commented  out  since 
the  default  output  file  is  LKCONBIF.DAT 

if  (argc  !=  2) 

printf ("Usage :  Prograjn  OutputFile  .. 

Now  exiting  system. .. \n" ; 
exit(l) 

if  ((outf ile=fopen(argv[l] , "a+t") )==NULL 
outf ile=fopen(argv[l] ,"at") 


/*  Input  default  values  for  constants  from  LKCONBIF.CFG  */ 

cfgfile  =  f openC'lkconbif . cfg" ,  "r"); 

/*  open  the  config  file  for  reading  */ 


if  (cfgfile  ==  NULL){ 

printf ("LKONBIF. CFG  file  not  found. .. .terminating  \n"); 
abort  0  ; 

} 


fscanf (cfgfile, 
fscanf (cfgfile , 
fscanf (cfgfile , 
fscanf (cfgfile , 
fscanf (cfgfile , 
fscanf (cfgfile , 
fscanf (cfgfile , 
fscanf (cfgfile , 


"If",  &BEGIN); 
"“/.f",  &END); 
"y.f",  &EXTIME); 
"%i",  &NP0INTS); 
"11",  &KMAX); 
"y.f",  &HTAMIN); 
"y.f",  &HTAMAX); 
"y.i",  &HTASTEP); 
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/*f  scanf  (cfgf  ile,  fetao)  ;  */ 


f scanf (cfgf ile ,  "%f", 
f scanf (cfgf ile ,  ""/of", 
f  scanf  (cfgf  ile ,  "yof", 
f  scanf  (cfgf  ile ,  "'/of", 
f  scanf  (cfgf  ile ,  "'/of", 
f  scanf  (cfgf  ile ,  "'/oe", 
f  scanf  (cfgf  ile ,  "'/e", 
f  scanf  (cfgf  ile ,  "'/oi", 
f  scanf  (cfgf  ile ,  "'/of", 
f  scanf  (cfgf  ile ,  "'/oi". 


&pi) : 

&alf a) ; 
&ext_length) ; 
fedelta) ; 
fephipcm) ; 
&photon_lif e) ; 
&carrier_life) ; 
femodel) ; 

&Y); 

feusernd) ; 


tao  =  carrier_life/photon_lif e; 

taoex  =  (2.0  *  ext_length)/(3.0el0  *  photon_lif e) ; 


printf("\n  BEGIN:  '/f", BEGIN  ); 
printf("\n  END:  '/of",END  ); 
printf("\n  EXTIME:  '/of", EXTIME  ); 
printf("\n  NPOINTS:  '/i",  NPOINTS) ; 
printf("\n  KMAX:  '/oi".KMAX  ); 
printf("\n  HTAMIN:  '/.f",  HTAMIN  ); 
printf("\n  HTAMAX:  If",  HTAMAX) ; 
printf("\n  HTASTEP:  '/i"  ,HTASTEP  ); 
printf("\n  tao:  '/of",tao  ); 
printf("\n  pi:  '/of",pi  ); 

printf("\n  alfa:  '/f",  alfa  ); 
printf("\n  taoex:  '/f", taoex  ); 
printf("\n  delta:  '/f ", delta) ; 
printf("\n  phipcm:  '/f",  phipcm) ; 
printf("\n  photon_life:  '/e",  photon_life) ; 
printf("\n  carrier_lif  e :  '/e",  carrier_lif  e) ; 
printf("\n  Y:  If",  Y) ; 
printf("\n  model:  '/i",  model); 
printf("\n  Random:  '/i  \n" ,  usernd)  ; 

if  (NPOINTS  ==  0){ 

printf  ("Values  not  properly  read  into  pr ogr aim.  ..  .terminating  \n"); 
abort  0 ; 

} 


outfile  =  f openC'lkconbif .dat" ,  "w"); 
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/* 


*/ 


1^ - i,! 

data.dim=4; 

if (model==0)  tau=taoex*sqrt (2 . O*pi/tao) ; 
if (model==l)  tau=taoex*sqrt(2.0*pi/tao) ; 
printf("\n  tau:  \n",  tau) ; 


dimen=4; 
taupi=tau; 
tpi=4.0*atan(l .0) ; 

/*  However  final  printed  time  is  in  units  of  tp  */ 

- */ 

h=(END-BEGIN)/(MPOINTS*KMAX) ; 
nstack=init_stack(tau,h, dimen) ; 
ystar=vector(l .data. dim) ; 
save=vector ( 1 , data . dim) ; 


/* 

/* 

/* 


- */ 

initial  conditions - */ 

- */ 


ystar[l]=0.10;  /*9.8;  .000645217;*/ 

ystar [2] =0 . 1 ;  /*  9.8;  0.0;  */ 

ystar[3]=0.01;  /*0.0;  -.014136;  */ 

ystar [4] =0.0; 

- -i,! 

/* - rkdumb - */ 

/* - */ 


f or (i=l ; i<=HTASTEP ; i++) { 

data. out=HTAMIM+(f loat) (i-1)* 

(HTAMAX- (HTAMIN) ) / ( (float) (HTASTEP-1) ) ; 
rkdumb (ystar , data . dim , BEGIN , END , NP0INTS*KMAX , der ivs ) ; 
for(j=l;  j<=data.dim; j++)  ystar [j]=save [j] ; 
if(usernd  ==  1){ 
if (model==0){ 

ystar [1] =6 . 0*(1 . 0-2 . 0*ranl (&iseed) ) ; 
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ystar  [2]  =2 . 0*tpi*raiil  (feiseed)  ; 
ystar [3] =6 . 0* ( 1 . 0-2 . 0*ranl (feiseed) ) ; 

} 

if (model==l){ 

ystar [1] =2 . 5*ranl (feiseed) ; 
ystar [2] =4 . 0*tpi*ranl (feiseed) ; 
ystar [3] =2 . 5*(1 . 0-2 . O+ranl (feiseed)) ; 

} 

/*printf  ("yl=%f  y2=y.f  y3=yof",  ystar[l],  ystar [2]  ,  ystar[3]);*/ 

} 

printf ("Iteration  #  completed  is  now  “/.d  \n",i); 

/*  if(i==l)  break;  */ 

} 


f closeCoutf ile) ; 
free_vector(save, 1 ,data.dim) ; 
free_ve ct or (ystar , 1 , data. dim) ; 
free_matrix(stack, 1 , dimen, 1 ,nstack+l) ; 
> 

#define  lA  16807 
#define  IM  2147483647 
#define  AM  (1.0/IM) 

#define  IQ  127773 

#define  IR  2836 

#define  NTAB  32 

#define  NDIV  (1+(IM-1) /NTAB) 

#define  EPS  1.2e-7 
#define  RNMX  (1.0-EPS) 

float  rani (long  *idum) 

{ 

int  j; 
long  k; 

static  long  iy=0; 
static  long  iv[NTAB]; 
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float  temp; 


if  (*idum  <=  0  I  I  ! iy)  { 

if  (-(*idum)  <  1)  *idum=l; 
else  *idum  =  -(*idum); 
for  (j=MTAB+7;j>=0;j— )  { 
k=(*idum)/IQ; 

*idum=IA*(*idum-k*IQ)-IR*k; 
if  (*idum  <  0)  *idum  +=  IM; 
if  (j  <  NTAB)  iv[j]  =  *idum; 

} 

iy=iv[0]  ; 

} 

k=(*idum)/IQ; 

*idum=IA*(*idum-k*IQ)-IR*k; 

if  (*idum  <  0)  *idum  +=  IM; 

j=iy/NDIV; 

iy=iv[j]  ; 

iv[j]  =  *idum; 

if  ((temp=AM*iy)  >  RNMX)  return  RNMX; 
else  return  temp; 

} 

#undef  lA 
#undef  IM 
#undef  AM 
#undef  IQ 
#undef  IR 
#undef  MTAB 
#undef  MDIV 
#undef  EPS 
#undef  RNMX 


void  rkdumb(vstart ,nvar ,xl ,x2 ,nstep,derivs) 
int  nvar.nstep; 
float  vstart [] ,xl ,x2; 
void  (*derivs)(); 
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{ 


int  i,k; 
float  x,h,s; 

float  *v,*vout,*dv,*vector() ; 
float  wplus , wminus , amp,laste, ampl ; 
void  rk4() ,nrerror() ,free_vector() ; 
v=vector(l ,nvar) ; 
vout=vector(l ,nvar) ; 
dv=vector ( 1 , nvar ) ; 
for  (i=l ; i<=nvar ; i++)  { 
v[i]=vstart  [i]  ; 

} 

x=xl ; 

h=(x2-xl)/nstep; 
for  (k=l ;k<=nstep ;k++)  { 
derivs(x, v,dv) ; 

rk4 ( V , dv , nvar , x , h , vout , der ivs ) ; 
if  (x+h  ==  x) 

nrerrorC'Step  size  too  small  in  routine  RKDUMB"); 

X  +=  h; 

s=(k+l)y.KMAX; 
for  (i=l ; i<=nvar ; i++)  -[ 
v[i]=vout  [i]  ; 

} 

/*  if(!s) 

fprintf  (outf  ile ,  "’/.f  ’/.f  %f  °/of  \n"  ,x,  v  [1]  ,v  [2]  ,  v  [3]  )  ;  */ 

} 

f or(i=l ; i<=nvar ; i++)  save[i]=v[i] ; 
free_vector(dv, 1 , nvar) ; 
free_vector(vout , 1 , nvar) ; 
free_vector(v, 1 ,nvar) ; 


/* - 

/* - rk4 


*/ 

*/ 

*/ 


void  rk4(y ,dydx,n,x,h,yout jderivs) 
float  y[]  ,dydx[]  ,x,h,yout  []  ; 
int  n; 

void  (*derivs)(); 

{ 
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int  i ; 

float  xh,hh,h6,*dym,*dyt ,*yt,*vector() ; 

dym=vector (1 ,n) ; 
dyt=vector (1 ,n) ; 
yt=vector(l ,n) ; 
hh=h*0.5; 
h6=h/6.0; 
xh=x+hh ; 

for  (i=l ; i<=n; i++)  yt [i] =y [i] +hh*dydx [i] ; 
derivs(xh,yt ,dyt) ; 

for  (i=l ; i<=n; i++)  yt  [i]=y[i]+hh*dyt [i] ; 
derivs (xh , yt , dym) ; 
for  (i=l ; i<=n; i++)  { 

yt [i] =y [i] +h*dym [i] ; 
dym[i]  +=  dyt  [i]  ; 

} 

derivs (x+h , yt , dyt ) ; 
for  (i=l ; i<=n; i++) 

yout  [i]  =y  [i]  +h6*  (dydx  [i]  +dyt  [i]  +2 . 0*dyin  [i] ) ; 
free_vector(yt , 1 ,n) ; 
free_vector(dyt , 1 ,n) ; 
free_vector(dym, 1 ,n) ; 


-*/ 

*/ 

-*/ 


/* - 

/* - derivs- 

jif. - 

void  derivs(t ,y,dy) 
float  t,y []  ,dy []  ; 

{ 

int  i ; 

float  hta,eps , fix, zZjtpi , gain; 
i=write_stack(y) ; 


hta=data. out ; 
eps=sqrt(tao/(2.0*pi)) ; 
tpi=4 . 0*atan(l . 0) ; 

gain=l . 0/sqrt ( 1 . 0+Y*pi* (y [1] *y  [1] -1.0)); 
fix=stack[2] [i]+y[2]-phipcm-2*delta*(t-taupi/2.0) ; 
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fix=fix-2.0*tpi*(float) (int) (fix/ (2.0*tpi)) ; 
if (model==l){ 

dy [1] =0 . 5* ( (2 . 0*y [3] +eps) *gain-eps) *y [1] + 
eps*hta*stack[l] [i]*cos(fix) ; 

dy[2]=  alfa*y[3]-eps*hta*(stack[l]  [i] /y  [1]  )*sin(f  ix) ; 

dy [3] =  (pi- (y [3] /eps) -pi* ( 1 . 0+ (2 . 0*y  [3] /eps ) ) * 
gain*y[l]*y[l])/(2.0*pi) ; 
zz=(y[l]-l  .0) ; 

} 

if (model==0){ 

dyCl]  =v[3] +hta*alf a*cos(f ix)/sqrt(2*(pi/tao) ) ; 
dy[2]=y[3]  ; 

dy [3]  =-y [1] -(l+2*pi)*y [3] /sqrt (2*pi*tao) ; 
dy[4]=y[2]  ; 
zz=y[l]  ; 

} 

if (model==0){ 

if  ( (t>=EXTIME)&&(lastde>=0 . 0)&:&(dy [1]  <=0 . 0) )  fprintf  (outf ile , 

"%f  y,f  y«f\n",  data.out*alfa/(sqrt(2.0*pi/tao)) ,  zz*alfa,  t) ; 
if ( (t>=EXTIME)&&(lastde<=0 . 0)&&(dy [1] >=0 . 0) )  fprintf (outf ile , 

"y,f  y.f  y.f \n"  ,  data.out*alfa/(sqrt(2.0*pi/tao)),  zz*alfa,  t) ; 
f f lush(outf ile) ; 
lastde=dy[l]  ; 

} 

if (model==l){ 

if ( (t>=EXTIME)&&(lastde>=0 . 0)&&(dy  [1] <=0 . 0) )  fprintf (outf ile , 

"y.f  "/of  yof\n",  data. out,  zz,  t)  ; 

if ( (t>=EXTIME) &&(lastde<=0 . 0)&&(dy [1] >=0 . 0) )  fprintf (outf ile , 

"y,f  "/of  y.f\n",  data-out,  zz,  t) ; 

ff lush(outf ile) ; 
lastde=dy[l] ; 

} 


/* - Outputs  phase - */ 

/*if ( (t>=EXTIME)&&(lastde>=0 . 0)&&(dy [2] <=0 . 0) ) 


119 


fprintf  (outf  ile,  "'/of  '/f  '/of\n", 

data. out*alf a/ (sqrt(2*pi/tao)) ,y[2]  ,t) ; 
if ( (t>=EXTIME)&&(lastde<=0 . 0)&&(dy [2] >=0 . 0) ) 
fprintf  (outf  ile,  "'/of  %f  '/f \n"  , 

data.out*alfa/(sqrt(2*pi/tao)) ,y[2] ,t) ; 
ff lush(oiitf  ile)  \ 
lastde=dy [2] ;  */ 

return; 

} 

/* - current - */ 


- */ 

int  init_stack(t , step, dim) 
float  tjStep; 
int  dim; 

{ 

int  i,j; 

nstack=(int) (4 . 0*t/step) ; 
stack=matrix(l ,dim, 1 ,nstack+l) ; 

f or(i=l ; i<=dim; i++) 

f or( j=l ; j<=nstack; j++)  stack[i] [j]=0.0; 
return  nstack; 

} 

int  write_stack(x) 
float  x[]  ; 

{ 

static  int  j=l; 
stack[l]  [j]=x[l]  ; 
stack  [2]  [j]=x[2]  ; 
if ( j==nstack)  j=0; 

j++: 

return  j ; 

} 

/* - */ 

void  nrerror(error_text) 
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char  error_text  [50] ; 

{ 

void  exitO  ; 

printf ("Numerical  recipes  run-time  error ... \n") ; 
printf  (""/osXn"  ,error_text) ; 
printf (".. .now  exiting  to  system  \n"); 
exit(l) ; 

} 

float  *vector(nl ,nh) 
int  nl,nh; 

/*  Allocates  a  float  vector  with  range  nl...nh  */ 
float  *v; 

v=(float  *)malloc( (unsigned)  (nh-nl+l)*sizeof (float) ) ; 
if(!v)  nrerrorC'allocation  failure  in  vectorO"); 
return  v-nl ; 

} 

void  free_vector(v,nl,nh) 
float  *v; 
int  nl ,  nh ; 

/*  Frees  a  float  vector  allocated  by  vectorO  */ 

{ 

free ((char*)  (v+nl)); 

} 

float  **matrix(nrl ,nrh,ncl ,nch) 
int  nrl ,nrh,ncl ,nch; 

{ 

int  i ; 
float  **m; 

/*  Allocate  pointers  to  rows  */ 

m=(float  **)malloc( (unsigned)  (nrh-nrl+l)*sizeof (float*) ) ; 
if(!m)  nrerrorC'Allocation  failure  1  in  matrix"); 
m-=nrl ; 

/*  Allocate  rows  and  set  pointers  to  them  */ 
f or (i=nrl ; i<=nrh; i++)  { 

m[i]=(float  *)malloc( (unsigned)  (nch-ncl+l)*sizeof (float)) ; 
if(!m[i])  nrerrorC'allocation  failure  2  in  matrix"); 
m[i]  -=ncl ; 
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} 

/*  Return  pointer  to  array  of  pointers  to  rows  */ 
return  m; 

} 

void  f ree_matrix (m , nrl , nrh , ncl , nch) 
float  **m; 

int  nrl, nrh, ncl, nch; 

/*  Frees  a  matrix  allocated  with  matrixO  */ 

{ 

int  i ; 

for(i=nrh;i>=nrl;i — )  free((char*)  (m[i]+ncl)); 
free((char*)  (m+nrl)); 

} 


Appendix  C.  LKC0N2.  C  Source  Code 

f* - lkcon2.  c - */ 

/* - Intergrates  the  conjugate  Lang  Kobayashi - */ 

/* — model  of  semiconductor  equations.  - */ 

/*  including  the  effects  of  gain  saturation  */ 

/* — This  is  a  code  for  delay-differential  equations - */ 

/* - with  the  appropriate  scalings  perfomed - */ 

f* - 10  Aug  95 - */ 

/* - Original  LICON.C  prograim  developed  by  */ 

/* - Dr  Gavrielides  Phillips  Laboratory  */ 

/*  modified  on  6  Dec  96  by  Gordon  T.  Hengst  */ 

/*  constants  moved  to  a  config  file  so  that  the  code  does  not  */ 
/*  have  to  be  re-compiled  for  every  constant  change .  */ 

/*  modified  on  22  April  97  to  include  the  effects  of  gain  */ 
/*  saturation  */ 


#include  <math.h> 
#include  <stdio.h> 
#include  <malloc.h> 


/*These  #define  statements  are 


replaced  with  variable  statements 


/*  #define 
/*  #define 
/*  #define 
/*  #define 
/*  #define 
/*  #define 
/*  #define 
/*  #define 


BEGIN  0.0  */ 

END  500.0  */  7*2000.0  */ 

EXTIME  400.0  */  7*985.0*7 

NPOINTS  32768  *7 
KM AX  1  *7 

HTAMIN  .00018  *7  7*. 0017 

HTAMAX  .001  *7  7*. 0006 

HTASTEP  100  *7 


*7 

*7 


*7 


7*  #define  inject  .0006425  *7 

7*  #define  tao  1428.57  *7 

7*  #define  pi  2.33  *7 

7*  #define  alfa  4.0  *7 

7*  #define  taoex  1095.23  *7  7*132.46*7 

7*  #define  delta  0.0  *7 

7*  +delta  is  injection  below  the  slave  frequency  *7 
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/*  #define  phipcm  0.0  */ 

/*  #define  model  1  */  /*  model=0  phase  equation  */ 

/*  model=l  full  model  */ 


float  BEGIN,  END,  EXTIME; 
int  NPOINTS,  KMAX; 
float  HTAMIN,  HTAMAX; 
int  HTASTEP; 

float  inject,  tao,  pi,  alfa,  taoex,  delta,  phipcm,  Y; 
int  model,  usernd; 

void  rkdumbO ; 
void  rk4(); 
float  *vector(); 
void  free_vector 0 ; 
void  nrerrorO; 
float  **matrix(); 
void  free_matrix() ; 
int  init_stack() ; 
int  write_stack() ; 
float  ranlO; 
void  f ourl 0 ; 


float  **y=0,*xx=0; 
float  **stack,*save; 

float  lastde=0 . 0 ,taupi , time2 ,photon_lif e , 
carrier_life ,ext_length,dir2; 
FILE  *outfile,*fourrf ,*fourfm; 

FILE  *cfgfile; 


int  nstack, iseed=-872581 , inti ; 
struct  given  { 

int  dim; 
float  out; 

}  data; 
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main(void) 


{ 

void  derivsO  ; 

void  exit  0 ; 

float  *ystar; 

register  int  i; 

float  h,tau, time, SCALE,  tpi; 

int  dimen, j,n; 


/„ - */ 

fourrf=f openC'f ourierrf .dat" , "wt") ; 
if (f ourrf==NULL)  { 

printfC'f ourierrf .dat  file  has  not  been  opened"); 
exit  (1); 

> 

fourfm=f openC'f ourierfm. dat" , "wt") ; 
if (f ourfm==NULL)  { 

printf ("fourierfm.dat  file  has  not  been  opened"); 
exit  (1) ; 

} 


/Ht  ***  Input  default  values  from  configuration  file  */ 

cfgfile  =  f openC'lkcon. cfg" , "r") ; 

/*  open  the  config  file  for  reading  */ 
if  (cfgfile  ==  NULL){ 

printf ("LKCON. CFG  file  not  found  ...  terminating  \n"); 
abort  0 ; 
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fscanf  (cfgfile,  "U"  ,  &BEGIN)  ; 
fscanf  (cfgfile,  &END)  ; 

fscanf  (cfgfile,  '"/.f",  &EXTIME) ; 
fscanf  (cfgfile,  ""/.i",  &MPOINTS)  ; 
fscanf  (cfgfile ,  "“/oi",  &KMAX) ; 
fscanf  (cfgfile,  ""/.f",  &HTAMIN)  ; 
fscanf  (cfgfile,  "y.f",  &HTAMAX) ; 
fscanf (cfgfile,  "%i",  &HTASTEP) ; 
fscanf (cfgfile,  "%f",  feinject); 

/*fscanf  (cfgfile,  "“/of",  &tao) ;  */ 

fscanf  (cfgfile,  "“/of",  &pi) ; 
fscanf  (cfgfile ,  "“/of",  &alf  a)  ; 
fscanf  (cfgfile ,  "“/of",  &ext_lengt]i)  ; 
fscanf  (cfgfile ,  "“/of",  fedelta)  ; 
fscanf  (cfgfile ,  "“/of",  fephipcm) ; 
fscanf  (cfgfile ,  "“/e",  &photon_lif  e)  ; 
fscanf  (cfgfile ,  "“/e",  &carrier_lif  e)  ; 
fscanf  (cfgfile ,  "“/oi",  femodel) ; 
fscanf  (cfgfile ,  "“/of",  &Y) ; 
fscanf  (cfgfile ,  "“/oi",  feusernd) ; 

tao  =  carrier_lif e/photon_lif e ; 

taoex  =  (2*ext_length)/(3el0  *  photon_life) ; 


printf("\n  BEGIN:  “/f", BEGIN  ); 
printf("\n  END:  ’/of", END  ); 
printf("\n  EXTIME:  “/f", EXTIME  ); 
printf("\n  NPOINTS:  ’/oi"  ,NP0INTS  ); 
printf("\n  KMAX:  “/oi",KMAX  ); 
printf("\n  HTAMIN:  “/of",  HTAMIN  ); 
printf("\n  HTAMAX:  “/of",  HTAMAX) ; 
printf("\n  HTASTEP:  “/i"  ,HTASTEP  ); 
printf("\n  inject:  “/of",  inject); 
printf("\n  tao:  “/of", tao  ); 
printf("\n  pi:  “/of", pi  ); 
printf("\n  alfa:  “/of",alfa  ); 
printf("\n  taoex:  “/of", taoex  ); 
printf("\n  delta:  ’/of", delta  ); 
printf("\n  phipcm:  “/of",  phipcm) ; 
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printf("\n  photon_life:  "/oe",  photon.lif  e)  ; 
printf("\n  carrier_lif e :  yoe",  carrier_life) ; 
printf("\n  Y:  y„f",  Y) ; 
printfC'Xn  model:  y„i",  model); 
printf("\n  Random:  "/.i  \n" ,  usernd) ; 

time2  =  1 . 0/photon_lif e ; 

if  (NPOINTS  ==  0)  { 

printf ("Values  not  properly  read  into  program  ...  terminating  \n"); 
abort  0 ; 

} 

outfile  =  fopenC'lkcon.dat" ,  "wt"); 

- */ 

data.dim=4; 

if (model==0)  tau=taoex*sqrt(2.0*pi/tao) ; 
if (model==l)  tau=taoex*sqrt(2.0*pi/tao) ; 

printfC'Xn  tau:  “/f  Xn" ,  tau) ; 

dimen=4; 

taupi=tau; 

tpi=4 . 0*atan(l . 0) ; 

/*  However  final  printed  time  is  in  units  of  tp  */ 

- 

h=(END-BEGIN)/(NPOINTS*KMAX) ; 
nstack=init_stack(tau,h, dimen) ; 
ystar=vector(l , data. dim) ; 
save=vector(l, data. dim) ; 
xx=vector(l ,NP0INTS+1) ; 
y=matrix(l , data. dim, 1 ,NP0INTS+1) ; 


- 

/* - initial 

- 

if (usernd  ==0){ 
ystar[l]=  0.52667; 
ystar[2]=  1.300658; 


- ^1 

conditions - */ 

- if.! 

/*  use  constant  values  */ 
7*9.8;  .000645217;*/ 

/*  9.8;  0.0;  */ 
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ystar[3]=  -0.265671;  /*0.0;  -.014136;  */ 

ystar [4] =0.0; 

} 

if(usernd  ==1)  { 

ystar [1]  =  2 . 5*ranl (feiseed) ; 

ystar [2]  =  4.0*tpi*ranl(&iseed) ; 
ystar[3]  =  2 . 5*(1 . 0-2 . 0*ranl (feiseed)) ; 

} 

/* - */ 

/* - rkdumb - */ 

data.out=inject; 

rkdumb (ystar , data . dim , BEGIN , END , NP0INTS*KMAX , derivs) ; 

/* - */ 

/* — Use  this  part  of  code  only  to  catch  isolated  solutions - */ 

/* 

f or(i=l ; i<=HTASTEP ; i++){ 

data. out=HTAMIN+(float) (i-1)* 

(HTAMAX-(HTAMIN))/((float)(HTASTEP-l)) ; 

inti=i ; 

data . out=HTAMIN ; 

ystar [1] =ranl (feiseed) ; 

ystar [2] =2 . 0*tpi*ranl (feiseed) ; 

ystar [3] =1 . 0-2 . 0*ranl (feiseed) ; 

printf  ("ysl=y,f  ys2=yof  ys3=yof  \n"  , 

ystar [1] , ystar [2] , ystar [3] ) ; 
rkdumb (ystar , data . dim , BEGIN , END , NP0INTS*KMAX , derivs ) ; 
save [2] =save [2] -2 . 0*tpi*(int) (save [2] /(2 . 0*tpi) ) ; 
for(j=l ; j<=data.dim; j++)  ystar [j] =save [j]  ; 
ystar [l]=ranl(&iseed) ; 
ystar [2] =2 . 0*tpi*ranl (feiseed) ; 
ystar [3] =1 . 0-2 . 0*ranl (&iseed) ; 
printf  ("Iteration  #  completed  is  now  yod\n",i); 
if(i==l)  break:; 

} 

*/ 

/* - */ 

if (model==0)  SCALE=sqrt (tao/ (2 . 0*pi) ) ; 
if (model==l)  SCALE=sqrt(tao/ (2.0*pi)) ; 


for(i=l;i<=NP0INTS;i++) 
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fprintf  (outf  ile,  "'/of  “/f  '/f  '/of\n", 

xx[i]*SCALE,  y[l]  [i]  ,y[2]  [i]  ,y[3]  [i] ) ; 


time=(xx [1001] -xxElOOO] )*SCALE/time2; 


/* -  Fourier  Transform  of  the  Amplitude - */ 

free_vector(xx, 1 ,NP0INTS+1) ; 
xx=vector(i ,2*NP0INTS+2) ; 
if(model==l)  { 

for  (i=l;i<=2*NP0INTS-l;i=i+2)  xx[i]=y[l] [(i+l)/2] ; 
f or (i=2 ; i<=2*NP0INTS ; i=i+2)  xx [i] =0.0; 

} 

if(model==0)  { 

for  (i=l;i<=2*NP0INTS-l;i=i+2)  xx [i] =1 . 0+y [1]  [(i+1) /2] ; 
f or (i=2 ; i<=2*NP0INTS ; i=i+2)  xx [i] =0.0; 

} 

fourl(xx,NP0INTS,l) ; 

f  or  (i=NPQINTS+l ;  i<=2*NP0INTS-l ;  i=i+2)  f  printf  (f  ourrf  ,  "'/.f  '/.f  \n"  , 

(float) (-(2*NP0INTS+l-i)/2)/( (float) (MP0INTS)*time) , 

(xx [i] *xx [i] +XX [i+1] *xx [i+1] ) / ( (float )NP0INTS) ) ; 
f  or(i=l ;  i<=MP0INTS-l ;  i=i+2)  fprintf(f  ourrf ,  ""/.f  "/f  \n"  , 

(float) ((i-l)/2)/((float) (NP0INTS)*time) , 

(xx [i] *xx [i] +XX [i+1] *xx [i+1] )/( (f loat)NP0INTS) ) ; 
fclose(f ourrf ) ; 

/+ - Fourrier  Transform  of  the  FM  spectrum - */ 

if(model==l)  { 

f or (i=l ; i<=2*NP0INTS-l ; i=i+2) 

xx[i]=y[l]  [(i+l)/2] *cos(y [2]  [(i+l)/2] ) ; 
for(i=2;i<=2*NP0INTS;i=i+2)  xx[i]=y[l] [i/2] *sin(y [2] [i/2] ) ; 

} 

if(model==0)  { 

f or(i=l ; i<=2*NP0INTS-l ; i=i+2) 

xx[i]  =  (1.0+y[l]  [(i+l)/2]  )*cos(y  [2]  [(i+l)/2]); 
for(i=2;i<=2*NP0INTS;i=i+2)  xx[i]=(l .0+y [1] [i/2] )*sin(y [2] [i/2]) 
} 

fourl(xx,NP0INTS,l) ; 
f or (i=MP0INTS+l ; i<=2*MP0INTS-l ; i=i+2) 

fprintf  (f  ourfm,  "'/of  '/of\n", 

(float) (-(2*NP0INTS+l-i)/2)/((float)(NPQIMTS)*time) , 
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(xx [i] *XX [i] +XX [i+1] *xx [i+1] ) / ( (float) NPOIMTS) ) ; 
f  or(i=l ;  i<=MPOINTS-l ;  i=i+2)  fprintf  (f  ourfm,  "7«f  °/of\n", 

(float) ((i-l)/2)/((float)(NP0INTS)*time), 

(xx[i] *xx[i]+xx[i+l] *xx [i+1] )/( (float) NPOINTS)) ; 
f close (f ourfm) ; 

f close(outf ile) ; 
free_vector(save , 1 , data. dim) ; 
free_vector(ystar,l, data. dim) ; 
free_matrix(stack, 1 ,dimen, 1 ,nstack+l) ; 

} 


#define  SWAP(a,b)  tempr=(a)  ;  (a)=(b)  ;  (b)=tempr 

void  f ourl (data.nn, isign) 
float  data[]  ; 
int  nn,isign; 

{ 

int  n,mmax,m, j ,istep,i; 
double  wtemp.wr ,wpr,wpi ,wi , theta; 
float  tempr, tempi; 
n=nn«l ; 

f or (i=l ; i<n; i+=2)  { 
if(j>i)  { 

SWAP(data[j]  ,data[i]  )  ; 

SWAP(data[j  +  l]  ,data[i+l] ) ; 

} 

m=n>>l ; 

while (m>=2  &&  j>m)  { 
j-=m; 
m»=l; 

} 

j+=m; 

} 

mmax=2 ; 

while (n>mmax)  { 

istep=2*mmax; 

theta=6 . 283185307 17959/ (isign+mmax) ; 
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wtemp=sin(0 . 5*theta) ; 
wpr=-2 . 0*wtemp*wtemp ; 
wpi=sin(theta) ; 
wr=l .0; 
wi=0 . 0 ; 

f or(in=l  ;in<mmax;m+=2)  { 

f or(i=m; i<=n; i+=istep)  { 
j=i+mmax; 

tempr=wr*data[j] -wi*data [j+1] 
tempi=wr*data[j  +  l]  +wi*data[j] 
data [ j ] =data [i] -tempr ; 
data[j+l]  =data[i+l]  -tempi ; 
data [i] +=tempr ; 
data [i+1] +=tempi ; 

} 

wr= (wtemp=wr ) *wpr-wi*wpi+wr ; 
wi=wi*wpr+wtemp*wpi+wi ; 

} 

mmax=istep; 

} 

} 


#define  lA  16807 
#define  IM  2147483647 
#define  AM  (1.0/IM) 

#define  IQ  127773 

#define  IR  2836 

#define  NTAB  32 

#define  MDIV  (1+(IM-1)/NTAB) 

#define  EPS  1.2e-7 

#define  RNMX  (1.0-EPS) 

float  rani (long  *idum) 

{ 

int  j; 
long  k; 

static  long  iy=0; 
static  long  iv[NTAB] ; 
float  temp; 
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if  (*idum  <=  0  I  I  ! iy)  { 

if  (-(*idum)  <  1)  *idum=l; 
else  *idum  =  -(*idum); 
for  (j=NTAB+7;j>=0;j— )  { 
k=(*idum)/IQ; 

*idum=IA* (*idum-k*IQ) -IR*k ; 
if  (*idum  <  0)  *iduin  +=  IM; 
if  (j  <  MTAB)  iv[j]  =  *idum; 

} 

iy=iv[0]  ; 

} 

k=(*idum)/IQ; 

*idum=IA*  ( *iduiii-k*IQ)  -IR*k ; 

if  (*iduin  <  0)  *idum  +=  IM; 

j=iy/MDIV; 

iy=iv[j]  ; 

iv[j]  =  *idum; 

if  ((temp=AM*iy)  >  RNMX)  return  RNMX; 
else  return  temp; 

} 


#undef  lA 
#undef  IM 
#undef  AM 
#undef  IQ 
#undef  IR 
#undef  NTAB 
#undef  MDIV 
#undef  EPS 
#undef  RNMX 


void  rkdumb(vstart ,nvar,xl ,x2,nstep, derive) 
int  nvar,nstep; 
float  vstart [] ,xl ,x2; 
void  (*derivs)(); 

{ 

int  i,k; 
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float  x,h,s; 

float  *v,*vout,*dv,*vector() ; 
float  wplus , wminus ,cLmp,laste,ampl ; 
void  rk4() ,nrerror() ,free_vector() ; 
v=vector(l ,nvar) ; 
vout=vector(l ,nvar) ; 
dv=vector (1 ,nvar) ; 
for  (i=l ; i<=nvar ; i++)  { 
v[i]=vstart  [i]  ; 

} 

x=xl ; 

h= (x2-xl ) /nstep ; 
for  (k=l ;k<=nstep ;k++)  { 
derivs(x, v,dv) ; 

rk4 ( V , dv , nvar , x , h , vout , der ivs ) ; 
if  (x+h  ==  x) 

nrerrorC'Step  size  too  small  in  routine  RKDUMB") 

X  +=  h; 

s=(k+l)y.KMAX; 
if(!s)  xx[l+(k/KMAX)]=x; 
for  (i=l ; i<=nvar ; i++)  { 
v[i]=vout  [i] ; 

if(!s)  y[i] [l+(k/KMAX)]=vout[i] ; 

} 

/*  if(!s)  fprintf  (outf  ile,  "%f  %f  %f  '/.f  \n"  , 

x,v[l]  ,v[2]  ,v[3]  )  ; 

if  (( !  s)&&(x>EXTIME))  fprintf  (outf  ile ,  "/od  '/ofXn", 
inti ,dir2) ; 
ff lush(outf ile) ;*/ 

} 

for(i=l ; i<=nvar ; i++)  save[i]=v[i] ; 
free_vector(dv, 1 ,nvar) ; 
free_vector(vout , 1 , nvar) ; 
free_vector(v, 1 , nvar) ; 


/He - 

/* - rk4 

/=, - 


*/ 

*/ 

*/ 


void  rk4(y,dydx,n,x,h,yout ,derivs) 
float  y  []  ,dydx[]  ,x,h,yout  []  ; 
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int  n; 

void  (*derivs)(); 

{ 

int  i ; 

float  xh,]ih,h6,*dym,*dyt,*yt,*vector()  ; 

dym=vector (1 ,n) ; 
dyt=vector(l ,n) ; 
yt=vector (1 ,n) ; 
hh=h*0.5; 
h6=h/6.0; 
xh=x+hh ; 

for  (i=l ; i<=n; i++)  yt [i]=y[i]+hh*dydx[i] ; 
der ivs  (xh ,  yt ,  dyt ) ; 

for  (i=l ; i<=n; i++)  yt [i]=y[i]+hh*dyt [i] ; 
derivs (xh , yt , dym) ; 
for  (i=l ; i<=n; i++)  { 

yt  [i]  =y  [i]  +h*dym  [i]  ; 
dyin[i]  +=  dyt  [i]  ; 

} 

derivs (x+h , yt , dyt ) ; 
for  (i=l ; i<=n; i++) 

yout [i] =y [i] +h6* (dydx [i] +dyt [i] +2 . 0*dym [i] ) ; 
free_vector(yt , 1 ,n) ; 
free_vector(dyt , 1 ,n) ; 
free_vector (dym, 1 ,n) ; 

} 

1^ - */ 

/* - derivs - */ 

1^ - */ 

void  derivs(t ,y ,dy) 
float  t,y[]  ,dy[]  ; 

{ 

int  i ; 

float  hta,eps,fix,zz,tpi,  gain; 
i=write_stack(y) ; 


hta=data. out ; 
tpi=4.0*atan(l .0) ; 
eps=sqrt(tao/ (2.0*pi)) ; 
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gain=l . 0/sqrt (1 . 0+Y*pi*(y [1] *y  [1] -1 . 0)) ; 


fix=stack[2] [i] +y [2] -phipcm-2*delta*(t-taupi/2 . 0) ; 
fix=fix-2.0*tpi*(float) (int) (fix/ (2.0*tpi)) ; 
if (model==l){ 

dy [1] =0 . 5* ( (2 . 0*y [3] +eps) *gain-eps) *y [1] + 
eps*hta*stack[l] [i] *cos(f ix) ; 

dy[2]=  alfa*y[3]-eps*hta*(stack[l]  [i]/y[l]  )*sin(fix) 

dy  [3] =  (pi- (y [3] /eps) -pi* ( 1 . 0+ (2 . 0*y [3] /eps) ) * 
gain*y [1] *y [1] )/(2 . 0*pi) ; 
zz=(y[l]-l  .0) ; 

} 

if (model==0){ 

f ix=stack[2] [i] +y [2] -phipcm-2*delta*(t-taupi/2 . 0) ; 
f ix=f ix-2 . 0*tpi* (float) (int) (fix/ (2 . 0*tpi) ) ; 
dy [1] =y [3] +hta*alfa*cos(f ix)/sqrt(2*pi/tao) ; 
dy  [2]  =y  [3]  ; 

dy [3] =-y [1] -(l+2*pi)*y [3] /sqrt(2*pi*tao) ; 
dy  [4]  =y  [2]  ; 
zz=y[l]  ; 

} 

dir2=dy[2]  ; 

/*if ( (t>=EXTIME)&&(lastde>=0 . 0)&&(dy [1] <=0 . 0) ) 
f  printf  (outf  ile ,  "’/of  ‘/f  ”/of  \n"  , 

data. out*alf a/ (sqrt(2*pi/tao)) ,zz,t) ; 
if ( (t>=EXTIME)&&(lastde<=0 . 0)&&(dy  [1] >=0 . 0) ) 
fprintf  (outf  ile,  "/of  ’/f  '/f \n"  , 

data. out*alf a/(sqrt(2*pi/tao)) ,zz,t) ; 
ff lush(outf ile) ;  */ 

lastde=dy[l] ; 


return; 

} 

/* - current - */ 

- 
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int  iiiit_stack(t, step, dim) 
float  t,step; 
int  dim; 

{ 

int  i,j; 

nstack=(int) (4 . 0*t/step) ; 
stack=matrix(l ,dim, 1 ,nstack+l) ; 
for(i=l ; i<=dim; i++) 

for(j=l; j<=nstack; j++)  stack[i] [j]=0.0; 
return  nstack; 

} 

int  write_stack(x) 
float  x[]  ; 

{ 

static  int  j=l; 
stack[l]  [j]=x[l]  ; 
stack [2]  [j]=x[2]  ; 
if ( j==nstack)  j=0; 

j++; 

return  j ; 

} 

jif. - if.1 

void  nrerror(error_text) 
char  error_text  [50] ; 

{ 

void  exit () ; 

printf ("Numerical  recipes  run-time  error ... \n") ; 
printf  ("“/osXn"  ,error_text) ; 
printf (".. .now  exiting  to  system  \n"); 
exit(l) ; 

} 

float  ^vector (nl ,nh) 
int  nl,nh; 

/*  Allocates  a  float  vector  with  reinge  nl...nh  */ 

{ 

float  *v; 

v=(float  *)malloc( (unsigned)  (nh-nl+l)*sizeof (float)) ; 
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if(!v)  nrerrorC "allocation  failure  in  vectorO"); 
return  v-nl; 

} 

void  free_vector(v,nl,nh) 
float  *v; 
int  nl,nh; 

/*  Frees  a  float  vector  allocated  by  vectorO  */ 

{ 

free ((char*)  (v+nl)); 

} 

float  **matrix(nrl,nrh,ncl,nch) 
int  nrl,nrh,ncl,nch; 

{ 

int  i ; 
float  **in; 

/*  Allocate  pointers  to  rows  */ 

m=(float  **)malloc( (unsigned)  (nrh-nrl+l)*sizeof (float*)) ; 
if(!m)  nrerrorC'Allocation  failure  1  in  matrix"); 
m-=nrl ; 

/*  Allocate  rows  and  set  pointers  to  them  */ 
for(i=nrl; i<=nrh; i++)  { 

m[i]=(float  *)malloc( (unsigned)  (nch-ncl+l)*sizeof (float)) ; 
if(!m[i])  nrerror ("allocation  failure  2  in  matrix"); 
m[i]  -=ncl ; 

} 

/*  Return  pointer  to  array  of  pointers  to  rows  */ 
return  m; 

} 

void  f ree_matrix(m,nrl ,nrh,ncl ,nch) 
float  **m; 

int  nrl ,nrh,ncl ,nch; 

/*  Frees  a  matrix  allocated  with  matrixO  */ 

{ 

int  i ; 

f or(i=nrh; i>=nrl ; i — )  free((char*)  (m[i]+ncl)); 
free((char*)  (m+nrl)); 

} 
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